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ABSTRACT 


An existing FORTRAN IV computer program for blade- to-blade aerodynamic anal- 
ysis of turbomachine blades was revised to obtain a simpler program consistent with re- 
lated programs. The analysis is for two-dimensional, subsonic, compressible (or in- 
compressible), nonviscous flow in a circular or straight infinite cascade of blades, 
which may be fixed or rotating. The flow may be axial, radial, or mixed, and the 
stream channel thickness may change in the through-flow direction. The results include 
streamline coordinates, velocity magnitude and direction throughout the passage, and 
the blade- surface velocities. This report includes a complete description of the input 
required by the program and the program listing. 



REVISED FORTRAN PROGRAM FOR CALCULATING VELOCITIES AND STREAMLINES 
ON A BLADE-TO-BLADE STREAM SURFACE OF A TURBOMACHINE 
by Theodore Katsanis and William D. McNally 
Lewis Research Center 

SUMMARY 

An existing FORTRAN IV computer program for blade- to-blade aerodynamic analysis 
of turbomachine blades was revised to obtain a simpler program consistent with related 
programs. The analysis is for two-dimensional, subsonic, compressible (or incompres- 
sible), nonviscous flow in a circular or straight infinite cascade of blades, which may be 
fixed or rotating. The flow may be axial, radial, or mixed, and the stream channel 
thickness may change in the through-flow direction. 

The program input consists of blade and stream channel geometry, total flow condi- 
tions, inlet and outlet flow angles, and blade-to-blade stream channel weight flow. The 
output includes blade-surface velocities, velocity magnitude and direction at all interior 
mesh points in the blade-to-blade passage, and streamline coordinates throughout the 
passage. 

This report includes a complete description of the input required by the program and 
the program listing. 


INTRODUCTION 

In the design of blade rows for compressors and turbines, it is desirable to obtain 
fluid velocities in the blade-to-blade passage and particularly on the blade surfaces. The 
trend to highly loaded blading results in widely spaced blades with less of the passage 
within a guided channel between blades. Stream filament techniques, applicable only 
within guided channels, can therefore no longer be used to obtain velocities over the en- 
tire blade surfaces. However, finite- difference methods can be used to obtain a solution 
of the stream- function differential equation in both the guided and unguided portions of the 
passage. 

Computer programs have been written which generate coefficients for the finite- 



difference equations, solve these equations, and differentiate the resulting values of 
stream function to obtain velocities throughout the blade- to- blade passage and on the blade 
surfaces. This was done in reference 1 for single blade row turbomachines or cascades, 
and in reference 2 for tandem or slotted blade machines. 

When the program of reference 2 (TANDEM) was written, many improvements were 
made over the single blade row program (2DCP) of reference 1. This report describes a 
new program (TURBLE) which solves the same problem as 2DCP but incorporates all of 
the improvements of TANDEM. The coding in TURBLE is both simpler and more fool- 
proof than that of 2DCP. Another reason for writing TURBLE has to do with the magnifi- 
cation program (MAGNFY) described in reference 3. The input to TURBLE has the same 
form as the input to TANDEM, and this simplifies the input coding to the magnification 
program (MAGNFY) of reference 3. It also allows the person using both TURBLE and 
TANDEM to put his input data in the same form in both cases. Further, TURBLE allows 
more interior mesh points in the solution region, and has its own error package indepen- 
dent of the Lewis computer system. Finally, the output of TURBLE has been expanded 
and clarified compared to the output of 2DCP. 

Like 2DCP, TURBLE obtains the numerical solution for ideal, subsonic, compres- 
sible (or incompressible) flow for an axial-, radial-, or mixed- flow cascade of turboma- 
chine blades. The cascade may be circular or straight (infinite), and may be fixed or ro- 
tating. The coordinates used are meridional streamline distance and angle in radians. 

This report includes a complete description of input and program listing for 
TURBLE. The mathematical analysis, the detailed program procedure, and the program 
output for TURBLE are all very similar to that for TANDEM (ref. 2). 

A TURBLE source deck on tape is available from COSMIC (Computer Software Man- 
agement and Information Center), Computer Center, University of Georgia, Athens, 
Georgia 30601. The program number is COSMIC number LEW- 10788. 


SYMBOLS 

m meridional streamline distance, meters, see figs. 1 and 2 
r radius from axis of rotation, meters 
s angular blade spacing or pitch, rad 

Vg tangential component of absolute fluid velocity, meters/sec 
W fluid velocity relative to blade, meters/sec 

z axial coordinate, meters 

a angle between meridional streamline and axis of rotation, rad, see fig. 1 
/3 angle between relative velocity vector and meridional plane, rad, see fig. 1 
2 



0 relative angular coordinate, rad, see fig. 1 

X prerotation , meters^/sec 

^ ^'in 
2 

p density, kg/meter 

(jj rotational speed, rad/sec 

Subscripts: 

cr critical velocity 

in inlet or upstream 

le leading edge 

out outlet or downstream 

te trailing edge 


DESCRIPTION OF INPUT AND OUTPUT 

The computer program requires as input a geometrical description in m- 6 coordi- 
nates of the blade surfaces, a description in m-r coordinates of the stream channel 
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Figure 2. - Blade-to -blade surface of revolution. 
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Figure 3. Input form. Card column numbers appear at top. 




through the blades, appropriate gas constants, and operating conditions such as inlet 
temperature and density, inlet and outlet flow angles, weight flow, and rotational speed. 
Figures 1 and 2 show the m- 6 coordinate system for a typical blade- to-blade surface of 
revolution. Output obtained from the program includes velocity magnitude and direction 
at all interior mesh points in the blade-to-blade passage, blade- surface velocities, 
stream- function values throughout the blade-to-blade region of solution, and streamline 
locations. 


Input 

Figure 3 shows the input variables as they are punched on the data cards. There 
are two types of variables, geometric and nongeometric. The geometric input variables 



Figure 4, - Geometric Input variables on blade-to-blade solution region. 



Figure 5. - Geometric input variables on a blade. BETI and BETO angles must be given as true angle 
not as angles measured in m-0 plane. Use tan p = r d0/dm to obtain p, or measure true angle. 
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Figure 6. - Geometric input variables describing stream channel in meridional plane. 


are shown in figures 4 to 6. All input variables are described in the section which fol- 
lows, Further explanation of key variables is given in the section Instructions for Pre- 
paring Input. 


The input variables are as follows: 


GAM 

AR 

TIP 

RHOIP 

WTFL 

OMEGA 


specific heat ratio 

gas constant, joule/(kg)(°K) 

inlet stagnation temperature, °K 

O 

inlet stagnation density, kg/meter 

mass flow per blade for the stream channel, kg/sec 

rotational speed, w, rad/sec (Note that lo is negative if rotation 
is in the opposite direction of that shown in fig. 1) 


ORF value of overrelaxation factor to be used in the solution of the inner 

iteration simultaneous equations (If ORF = 0, the program 
calculates an estimated value for the overrelaxation factor. 

See p. 10 for discussion. ) 


BETAI 

BETAO 

CHORDF 

STGRF 

MBI 


inlet flow angle /3, along BG with respect to m-direction, deg 
(fig. 4) 

outlet flow angle along CF with respect to m-direction, deg 
(fig. 4) 

overall length of blade in m-direction, meters (fig. 4) 

angular 0- coordinate for center of trailing- edge circle of blade 
with respect to the center of leading- edge circle of blade, rad 
(fig. 4) 

number of vertical mesh lines from AH to BG inclusive (fig. 4) 


6 



MBO 


MM 

NBBI 

NBL 

NRSP 

RI1,RI2 

R01,R02 

BETI1,BETI2 


BET01,BET02 


SPLN01,SPLN02 


MSP1,MSP2 


THSP1,THSP2 


MR 


number of vertical mesh lines from AH to CF inclusive (fig. 4) 

total number of vertical mesh lines in the m-direction from AH to 
DE, maximum of 100 (fig. 4) 

number of mesh spaces in 0-direction between AB and GH, maxi- 
mum of 50 (fig. 4) 

number of blades 

number of spline points for stream channel radius (RMSP) and 
thickness (BESP) coordinates, maximum of 50 (fig. 6) 

leading- edge radii of the two blade surfaces, meters (fig. 5) 

trailing- edge radii of the two blade surfaces, meters (fig. 5) 

angles (with respect to m-direction) at tangent points of leading- 
edge radii with the two blade surfaces, deg (fig. 5) (These must 
be true angles in degrees. If angles are measured in the m- 0 
plane, i. e. d0/dm, BETH and BETI2 can be obtained from the 
relation tan ^ = r(d0/dm). ) 

angles (with respect to m-direction) at tangent points of trailing- 
edge radii with the two blade surfaces, deg (fig. 5) (These must 
also be true angles in degrees, like BETH and BETI2. ) 

number of blade spline points given for each surface as input, 

maximum of 50 (These include the first and last points (dummies) 
that are tangent to the leading- and trailing- edge radii (fig. 5). ) 

arrays of m- coordinates of spline points on the two blade surfaces, 
measured from the blade leading edge, meters (fig. 5) (The first 
and last points in each of these arrays can be blank or have a 
dummy value, since these points are calculated by the program. 

If blanks are used, and the last point is on a new card, a blank 
card must be used. ) 

arrays of 0- coordinates of spline points corresponding to MSPl and 
MSP2, rad (fig. 5) (Dummy values are also used here in posi- 
tions corresponding to those in MSPl and MSP2. ) 

array of m- coordinates of spline points for the stream channel 
radii and the stream channel thicknesses, meters (fig. 6) (MR is 
measured from the leading edge of the blade. These coordinates 
should cover the entire distance from AH to DE, and may extend 
beyond these bounds. The total number of points is NRSP. ) 
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RMSP 


array of r- coordinates of spline points for the stream channel 
radii, corresponding to the MR array, meters (fig. 6) 

array of stream channel normal thicknesses corresponding to the 
MR and RMSP arrays, meters (fig. 6) 

The remaining variables, starting with BLDAT, are used to indicate what output is de- 
sired. A value of zero for any of these variables will cause the output associated with 
that variable to be omitted. A value of 1 will cause the corresponding output to be printed 
for the final outer iteration only; 2, for the first and final iterations; and 3, for all outer 
iterations. Care should be used not to call for more output than is really useful. The 
following list gives the output associated with each of these variables. 


BLDAT 


AANDK 

ERSOR 

STRFN 

SLCRD 

INTVL 

SURVL 


all geometrical information which does not change from iteration to iteration 
(i. e. , coordinates and first and second derivatives of all blade surface spline 
points; blade coordinates and blade slopes where vertical mesh lines meet 
each blade surface; radii and stream channel thicknesses corresponding to 
each vertical mesh line; m- coordinate, stream channel radius and thick- 
ness, and blade surface angles and slopes where horizontal mesh lines inter- 
sect each blade; and ITV and IV arrays (internal variables describing the lo- 
cation of the blade surfaces with respect to the finite difference grid). ) 

the coefficient array, the constant vector, and the indexes of all adjacent points 
for each point in the finite-difference mesh (This information is needed for 
debugging the program only. ) 

the maximum change in the stream function at any point for each iteration of the 
SOR equation, eq. (A 8), ref. 2 

value of the stream function at each unknown mesh point in the region 
streamline 6- coordinates at each vertical mesh line, and streamline plot 
velocity and flow angle at each iterior mesh point 

m-coordinate, surface velocity, flow angle, distance along surface, and 
W/W^^ based on meridional velocity components where each vertical mesh 
line meets each blade surface; m-coordinate, surface velocity, flow angle, 
distance along surface, and W/W^^ based on tangential velocity components 
where each horizontal mesh line meets each blade surface; and plot of blade- 
surface velocities against meridional streamline distance, meters. (It is sug- 
gested that SURVL=3 be used. This will give surface velocities after each 
outer iteration, so that satisfactory velocities may be obtained even when 
final convergence is not reached. ) 
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Instructions for Preparing Input 


Units of measurement. - The International System of Units (ref. 4) is used through- 
out this report. However, the program does not use any constants which depend on the 
system of units being used. Therefore, any consistent set of units may be used in pre- 
paring input for the program. For example, if force, length, temperature, and time are 
chosen independently, mass units are obtained from force = mass x acceleration. The 
gas constant R must then have the units of force times length divided by mass times 
temperature (energy per unit mass per degree temperature). Density is mass per unit 
volume, and weight flow is mass per unit time. Output then gives velocity in the chosen 
units of length per unit time. Since any consistent set of units can be employed, the out- 
put is not labeled with any units. 

Blade and stream channel geometry . - The upper and lower surfaces of the blade are 
each defined by specifying three things: leading- and trailing- edge radii, angles at which 
these radii are tangent to the blade surfaces, and m- and 0- coordinates of several points 
along each surface. These angles and coordinates are used to define a cubic spline curve 
fit (ref. 5) to the surface. The standard sign convention is used for angles, as indicated 
in figure 5. 

A cubic spline curve is a piecewise cubic polynomial which expresses mathematically 
the shape taken by an idealized spline passing through the given points. Reference 5 de- 
scribes a method for determining the equation of the spline curve. Using this method, 
few points are required to specify most blade shapes accurately, usually no more than 
five or six, in addition to the two end points. As a guide, enough points should be speci- 
fied so that a physical spline passing through these points would accurately follow the 
blade shape. This means that the spline points should be closer where there is large 
curvature and farther apart where there is small curvature. 

The coordinates for either surface of the blade are given with respect to the leading 
edge, with the leading edge of the blade being defined as the furthest point upstream. 

The mean stream surface of revolution (as seen in the meridional plane, fig. 6) and 
the stream channel thickness are also fitted with cubic spline curves. The m- coordinates 
for the mean stream surface are independent of the m- coordinates for blade surfaces. 

Inlet and outlet flow angles . - The values of and are given as average val- 
ues on BG and CF, respectively. If the flow is axial these flow angles are the same as 
the flow angles at AH and DE. If flow is radial or mixed, and these angles are not known 
on BG and CF, /3j^ and must be calculated by equation (B15) of reference 1 or 
equation (B14) of reference 2. 

Defining the mesh. - A finite- difference mesh is used for the solution of the basic 
differential equation. A typical mesh pattern is shown in figure 7. The mesh spacing and 
the extent of the upstream and downstream regions are determined by the values of MBI, 
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Figure 7. - Typicai mesh in blade-to-blade solution region. 


MBO, and MM of the input. The mesh spacing must be chosen so that there are not more 
than 2500 unknown mesh points. 

Values of MBI, MBO, and MM should be determined so that the mesh which results 
has blocks which are approximately square. To achieve this, a value for NBBI is first 
chosen arbitrarily (15 to 20 is typical). NBBI is the number of mesh spaces spanning the 
blade pitch, s, where s = 27 t/NBL. Dividing s by NBBI gives the mesh spacing, HT, in 
the e- direction in radians. Muliplying HT by an average radius (RMSP) of the stream 
channel gives an average value for the actual mesh spacing in the 6»- direction. The value 
of CHORD should then be used with this tangential mesh spacing to calculate the approxi- 
mate number of mesh spaces along the blade in the m-direction. This will give MBO 
once MBI is chosen. Generally, MBI is given a value of 10. MM, likewise, is usually 
given a value 10 more than MBO. 

Overrelaxation factor . - ORF is the over relaxation factor used in each inner itera- 
tion in the solution of the simultaneous finite difference equations. (See ref. 2, p. 101). 
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ORF may be set to zero, or some value between 1 and 2. ORF is usually given as zero 
for the initial run of a given blade geometry and mesh spacing (MBI, NBBI, etc. ). In 
this case the program uses extra time and calculates an optimum value for ORF. It does 
this by means of an iterative process, and on each iteration the current estimate of the 
optimum value for ORF is printed. The final estimate is the one used by the program for 
ORF. If the user does not change the mesh indexes MBI, MBO, MM, and NBBI between 
runs, even though blade geometry or other input does change, he may use this final es- 
timate of ORF in the input, saving the time used in its computation. In all cases, if ORF 
is not zero, it should have a value greater than 1 and less than 2. 

Actually, the value of ORF is not as critical as the user might think. It gets more 
critical as the optimum value gets close to 2. For any run of a given set of data, only 
small changes will occur in the rate of convergence in SOR as long as the difference 
2.0- ORF is within 10 percent of its optimum value. 

Format for input data. - All the numbers on the card beginning with MBI and on the 
card beginning with BLDAT are integers (no decimal point) in a 5-column field (see 
fig. 3). These must all be right adjusted. The input variables on all other data cards are 
real numbers (punch decimal point) in a 10- column field. 

Incompressible flow. - While the program is written for compressible flow, it can 
be easily used for incompressible flow. To do so specify GAM =1.5, AR = 1000, and 
TIP = 10® as input. This results in a single outer iteration of the program to obtain the 
stream function solution. 

Straight infinite cascade . - The program is as easily applied to straight infinite cas- 
cades as circular cascades. Since the radius and number of blades (NBL) for such a 
cascade would actually be infinite, an artificial convention must be adopted. The user 
should pick a value for NBL, for instance 20 or 30. Then, since the blade pitch sr 
(fig. 4) is known, an artificial radius can be computed from 

^ ^ NBL*(sr) 

277 

This r should be used to compute the 0- coordinates requires as input (THSPl, THSP2, 
and STGR) by dividing coordinates in the tangential direction by r. 

Axial flow. - For a two-dimensional cascade with constant stream channel thickness, 
constant values should be given for the MR, RMSP, and BESP arrays. Only two points 
are required for each of these arrays in this case. The two values of MR should be 
chosen so that they are further upstream and downstream than the boundaries AH and DE. 
The two values of RMSP and BESP should equal the constants r and b. 
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TABLE I. - SAMPLE OUTPUT 


AXIAL STATOR 

- MEAN SECTION 

- 

COMPRESSIBLE 





GAM 

AR 


UP 

RHOI P 

WTFL 

WTFLSP 


1.4000000 

287.05300 


288.15000 

1.225COOO 

0.3146000 

-0 

-0 

BETAI 

BET AO 


CHOROF 

STGRF 




0 

-67. COOOOO 


0.4265000E-01 

-0.1116150 




MBl MBO 

MM NBBI 

NBL 

NRSP 





15 32 -0 

0 47 20 

50 

2 






BLADE surface 

1 — UPPER SURFACE 





RI 1 

ROl 

8ETI 1 

BETOl 

SPLNOl 



0.3810000t-02 

0.8890000E-03 

28.300000 

-72.400000 

7-OOOCOOO 



MSPl ARRAY 







-0 

0.8575000E-02 

0. 1715000E-01 

0.2572500E-01 

0.343C000E-01 

0. 3858800E-01 

-0 

THSPl ARRAY 







-0 

0. 1 769000E-01 

0. 1538000E-01 

-0.5310000E-02 

-0.4654000E-01 

-0. 7400000E-01 

-0 


BLADE SURFACE 

2 — LOWER SURFACE 




RI2 

RD2 

BETI2 

6ETC2 

SPLN02 


0.3810000E-02 

0.8890000E-03 

-14.200000 

-56. 100000 

6.0000000 


MSP2 ARRAY 






-0 

0. 0575OOOE-O2 

0. 1715000L-01 

0.2572500E-01 

0.3430000E-01 

-0 

THSP2 ARRAY 






-0 

-0. 1562000E-01 

-0.2854000E-01 

-0.5070000E-01 

-O.025COOOE-O1 

-0 


MR ARRAY 

-0.5000000E-01 0-1000000 

RMSP ARRAY 

0.3302000 0-3302000 

8ESP ARRAY 

0.1016000 0.1016000 



BLDAT 

AANOK 

ERSOR 

STRFN 

SLCRD 

INTVL 

SURVL 

L 

1 

2 

2 

2 

2 

2 

3 


BLADE DATA AT INPUT SPLINE POINTS 


M 

0.20037E-02 
0.85750E-02 
0.17150E-01 
0.25725E-01 
0.34300E-01 
J 0.3858bE-01 
^ 0.42608E-01 


BLADE 


SURFACE 


I 


THETA 

0.10159E-01 

0.176908-01 

0.15380E-01 

-0.53100E-02 

-0.465A0E-01 

-0.7<r000E-01 

-0.11080 


DERIVATIVE 

1.63066 

0.60358 

-1.23217 

-3.72155 

-5.54459 

-7.92944 

-9.54694 


2ND DERIV 
-129.936 
- 182.661 
-245.503 
-335.112 
-90.0848 
-1022.25 
217.603 


M 

0.28754E-02 

0.85750E-02 

0.17150E-01 

0.25725E-01 

0.34300E-01 

0.41023E-01 

V. 


BLADE SURFACE 2 


theta 

0.11448 

0.11004 

0.97 124E-01 

0.74964E-01 

0.43164E-01 

0.12547E-01 


DERIVATIVE 

-0.76632 

-0.96494 

-2.04512 

-3.12744 

-4.32324 

-4.50684 


2ND DERIV. 

57.4409 

-127.138 

-124.799 

-127.637 

-151.265 

96.6479 



3 S LEADING EDGE B-G 
TRAILING EDGE C-F 


FREESTREAM 
VELOC I TY 
61*9304 
180.067 


MAXIMUM VALUE 
FOR RHO*H 
241.239 
241-239 


CRITICAL 

VELOCITY 

310.645 

310.645 


BOUNDARY A-H 
BOUNDARY D-E 


4 < 


CALCULATED PROGRAM CONSTANTS 


PITCH 

0.1256637 


ITMIN 

-17 


lambda 


HT 


HMl 


0.6283185E-02 0. 2508824E-02 


UMAX 

19 


NUMBER OF INTERIOR MESH POINTS = 


SURFACE BOUNDARY VALUES 
BV 


SURFACE 

1 

2 


0 . 

1.00000 




BLADE DATA AT INTERSECTIONS 

, OF vertical 

MESH LINES WITH 

BLADES 

BLADE SURFACE 1 


BLADE SURFACE 2 

M TV 

OTOMV 


TV 

OTOMV 

0 0 

O.IOOOOE 

11 

0. 12566 

-O.IOOOOE 

0.25088E-02 0.109666-01 

1.56401 


0.11482 

-1.10043 

0.50176E-02 0.14447E-01 

1.20260 


0.11291 

-0.71758 

0.75265E-02 0.16958E-01 

0.79070 


0.11099 

-0.84944 

0.10035E-01 0.U373E-01 

0.32903 


0.10850 

-1.15031 

STREAM SHEET COORDINATES AND 

THICKNESS TABLE 


IM M R 

SAL 


B 

DB/I 

1 -0.351246-01 0.33020 

-0 


0.10160 

-0 

2 -0.32615E-01 0.33020 

-0 


0.10160 

-0 

3 -0.30106E-01 0.33020 

-0 


0.10160 

-0 

4 -0.27597E-01 0.33020 

-0 


0. 10160 

-0 

5 -0.25088E-01 0.33020 

-0 


0.10160 

-0 


2 < 


IM 


IV ARRAY 


I 

21 

41 

61 

81 


ITV ARRAY 
1 2 


19 

19 

19 

19 

19 


BETA CORRECTED 
TO BOUNDARY 
0 

-67.0000 


BLADE 

SURFACE 

NO. 


TABLE I. -Continued. SAMPLE OUTPUT 


"m coordimates of intersections of horizontal mesh lines with blade 


MH array - BLADE SURFACE I 


0.6132E-03 
0,3582E-02 
0, 19066-01 
0.22016-01 
0.2A20E-01 

- BLADE SURFACE 2 

MH 

0.61326-03 
0. A7596-02 
0.U40E-01 
0. 1539E-01 
0.1850E-01 


RMH 

0.3302 

0.3302 

0.3302 

0.3302 

0.3302 

0.3302 


RMH 

0-3302 

0-3302 

0.3302 

0.3302 

0.3302 


BEH 

0. 1016 
0. 1016 
0. 1016 
0. 1016 
0. 1016 
0. 1016 


BEH 

0. 1016 
0. 1016 
0. 1016 
0- 1016 
0. 1016 


8ETAH 
90.000 
57.U4L 
25-052 
-29.610 
-40.095 
-46. 766 


BETAH 

-57.041 

-13.294 

-23.596 

-31.073 

-36.169 


OTOMH 

-4.6708 

-0.7156 

-1-3228 

-1.8250 

-2.2140 


THETA COOROMATES OF HORIZONTAL MESH LINES 

IT THETA 

-17 -0.10681 

-16 -0. 10053 

-15 -0.94248E-01 

-14 -0.87965E-01 

-13 -0.81681E-01 


l.OOOCO 
1.00000 
1.00000 
1. 00000 
l.OOOQO 




0-22687771 0.2769^251 0.32706^63 

0.72817850 0.77811576 0.82799505 


0.22687770 0.2769-V254 0.32706478 

0.72817869 0.77811572 0.82799502 


0.22678433 0-27685902 0.32699915 

0-72827139 0.77819968 0.82806225 


STREAM FN. THETA 

0- 4000000 0.46832566*01 

1.0000000 0.1222286 

0.4000000 0.468 325 5E-01 

1- 0000000 0.1222287 

0.4000000 0.468361 8E-0 I 

1. 0000000 0.1222321 


0.37723198 0.42742761 0.47763271 

0.87782843 0.92763250 0.97742604 


0-37723205 0.42742773 0-47763292 

0.87782862 0.92763270 0-97742602 


0.37719061 0.42741434 0.47764858 

0-87787237 0.92764822 0.97741155 


STREAM FN. THETA 

0.6000000 0.7187305E-01 
0.2000000 0.2175627E-01 
0.6000000 O.71073OOE-OI 
0.2000000 0.2175627E-01 
0.6000000 0.7186363E-01 
0.2000000 0.2176819E*0l 


m 
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0201 



0401 


o 
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VELnCITIES AT INTERIOR MESH POINTS 


1 


13 ^ 


VELOCITY 

61.703 

62-060 

62.173 

61-819 

ANGLEIOEG) 

-0.01 

-0.03 

O-Ol 

0-03 

VELOCITY 

61.750 
62.127 
62.125 

61.751 

ANGLEIOEG) 

-0.02 

-0.02 

0.02 

0.02 

VELOCITY 

61*818 

62.175 

62.059 

61.702 

VELOCITY 

61.703 

62.060 

62-173 

61.819 

ANGLEIOEG ) 
0.02 
0.04 
-0.02 
-n.04 

VELOCITY 

61.750 
62.127 
62. 125 

61.751 

ANGLEIOEG) 

0-03 

0.03 

-0.03 

-0.03 

VELOCITY 

61.816 

62.175 

62.059 

61.703 


ANGLEIOEG) 

-0.03 

-0.01 

0.03 

0-01 

VELOCITY 

61.898 

62.199 

61.981 

61.676 

ANGLEIOEG) 

-0.03 

-0.00 

0.03 

o.oi 

VELOCITY 

61.981 

62.199 

61.898 

61.678 

ANGLEIOEG) 

-0.03 

0.01 

0.03 

-0.01 

ANGLEIOEG) 

0.04 

0.02 

-0.04 

-0.02 

VELOCITY 

61.898 

62.200 

61.980 

61.676 

ANGLE lOEG) 
0.04 
0.01 
-0.04 
-0.01 

VELOCITY 

61.981 

62.199 

61.898 

61.678 

ANGLEIOEG) 

0-04 

-0.01 

-0.04 

0.01 


14 ITERATION NO. 


MAXIMUM RELATIVE CHANGE IN DENSITY = 0.577-V 


15 < 


M 

O 

0.2509E^02 
0.5018E-02 
0.7526E-02 
0. lOO^E-Ol 


SURFACE VELOCITIES BASED ON MERIDIONAL COMPONENTS 


VELOCITY 

0 

99.866 

114-60 

134-41 

158-08 


BLADE 

ANGLE(DEG) 

90-00 

27.31 

21.66 

14.63 

6-20 


SURFACE I 
SURF- LENGTH 
0 

0.4405E-02 
0. 7165E-02 
0.9807E-02 
0. 1236E-01 


W/WCR 

0 

0-3215 

0.3689 

0.4327 

0-5089 


VELOCITY 

0 

74.794 

63.478 

59.630 

58.400 


BLADE 

ANGLEIOEG) 

-90.00 

-19.97 

-13.33 

-15.67 

-20.80 


SURFACE 2 
SURF. LENGTH 
0 

0.4372E-02 
0.6959E-02 
0.9547E-02 
0. 1219E-01 


W/WCR 

0 

0.2408 

0.2043 

0.1920 

0.1080 


♦ 

♦ 

♦ 

♦ 

* 

* 

♦ 

♦ 


SURFACE VELOCITIES BASED ON TANGENTIAL 


COMPONENTS 


0.6132E-03 
0. 3582E-02 
0. 1906E-01 
0.2201E-01 


BLADE 

VELOCITY 

11.357 

76.923 

100.30 

240.79 

240.01 


SURFACE I 
ANGLEIOEGI 
90.00 
57.04 
25-05 
-29.61 
-40. 10 


W/WCR 

0.3656E-01 
0.2476 
0.3229 
0.7751 
0. 7726 


15 


M 

0.6132E-03 
0.4759E-02 
0. 1 140E-01 
0. 1539E-01 


BLADE 
VELOC ITY 
47.792 
48-156 
57. 949 
62.953 


surface 2 

ANGLE! DEG ) 
-57.04 
-13.29 
-23.60 
-31.07 


W/WCR 

0. 1538 
0. 1550 
0-1865 
0.2027 



TABLE!. - Concluded. SAMPLE OUTPUT 


BLADE SURFACE VELOCITIES 


15 K 


0, 

50. 100. 

150. 

200. 

250. 

300. 

0. !-♦ 

X 

0*1 

I 

X 

X 

X 

X 

1 I 

X 

1 

1 

X 

X 

1 X ♦ 

X 

X 

X 

X 

X 

1 X 

X 

X 

X 

1 

X 

IX X ♦ 

1 

X 

X 

X 

X 

X X 

X 

1 

1 

X 

X 

X X 

1 

1 

X 

1 

X 

X X X 

♦ X 

X 

L 

X 

X 

1 X 

1 

1 

1 

X 

O.OXOOl 


1-+ — 

1 

X 

X — 

X 

X X 

X 

X 

1 

X 

1 

X X 

1 

X 

X 

1 

X 

X X 1 

1 

♦ 1 

1 

1 

X 

1 I 

X 

X 

X 

X 

X 

X 0 1 

X 

!♦ 

X 

X 

X 

X X 

X 

X 

X 

X 

1 

X 1 

X 

1 

X 

X 

X 

X X X 

X 

X ♦ 

X 

X 

X 

X 0 X 

X 

X 

1 

X 

0.0200X 


1 

X 

-+--1 

X — 

X 

X 0 1 

1 

1 

X 

X 

X 

X X 

1 

X 

» 1 

1 

X 

X 0 X 

X 

X 

+ X 

X 

X 

X X 

X 

1 

* X 

X 

X 

X XX 

X 

1 

♦ 1 

X 

X 

X 0 1 

X 

X * 

X 

X 

X 

X 1 

X 

X 

X 

1 

X 

X ox 

1 

X ♦ 

1 

X 

X 

X xo 

1 

1 ♦ 

X 

1 

U . U 5UUI i X” A i 

X 

X 1 0 

1 

X* 

1 

1 

X 

X X 

1 

X 

1 

X 

X 

1 X 

X X 

X* 

1 

1 

X 

X 1 

0 1 

X* 

X 

X 

X 

X X 

X X 

X ♦ 

I 

X 

X 

1 1 

ox 

X ♦ 

X 

X 

1 

X X 

10 

X ♦ 

X 

X 

1 

X X 

X 0 

X ♦ 

X 

X 

1 

I X 

1 

1 ♦ 

X 

X 

0.04001 

1 1 — 

X 

-0--4X* 

1 

X-; 

I 

1 1 

X 

♦ 1 

X 

X 

X 

X X 

X 

♦ X 

X 

X 

1 

X X 

X 

X 

X 

1 

X 

X X 

X 

X 

X 

X 

X 

X 1 

1 

X 

X 

X 

X 

X 1 

X 

X 

X 

X 

X 

X X 

1 

X 

X 

X 

X 

X X 

1 

X 

1 

1 

X 

X X 

1 

X 

X 

X 

C. 0500X 

0. 

^ 1 ^ 1 
50. 100. 

— x---— 

150. 

200. 

250- 

300 


350- 

— I — 
I 

I 

1 

I 

X 

I 

I 

1 

i 

1 

1 

1 

1 

I 

1 

1 

1 

L 

I 

I 

I 

I 

I 

1 

I 

X 

X 

— X — 
I 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

1 

X 

350- 


400. 

I 

I 

X 

X 

I 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

I 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

400- 


450. 

X 

X 

X 

X 

1 

X 

X 

X 

X 

X 

X 

1 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

1 

X 

X 

X 

X 

X 

X 

X 

X 

X 

450. 


500. 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

X 

500. 


VELOCITYIW) vs. MERIDIONAL STREAMLINE DISTANCE(M) DOWN THE PAGE 


CO 


- BLADE SURFACE X, 
’ BLADE SURFACE X, 

- BLADE SURFACE 2, 

- BLADE SURFACE 2, 


BASED ON MERIDIONAL COMPONENT 
BASED ON TANGENTIAL COMPONENT 
BASED ON MERIDIONAL COMPONENT 
BASED ON TANGENTIAL COMPONENT 


16 TIME 


L.0575 min 



Output 


Sample output is given in table I for the axial-flow stator example of reference 1. 

The blade shape is shown in figure 7. Since the complete output would be lengthy, only 
the first few lines of each section of output are reproduced herein. Most of the output is 
optional, and is controlled by the final input card, as already described. In some in- 
stances output labels are simply internal variable names. 

Each section of the sample output in table I has been numbered to correspond to the 
following description: 

(1) The first output is a listing of the input data. All items are labeled as on the in- 
put form (fig. 3). 

(2) This is the output corresponding to BLDAT. (See the list of input variables. 

See ref. 2 for meaning of undefined labels. ) 

(3) The relative free- stream velocity W; the relative critical velocity W ; and the 
maximum value of the mass flow parameter pW (corresponding to W = W ) are given at 
the leading edge of the blade (BG) and the trailing edge of the blade (CF). The inlet (out- 
let) free-stream flow angle /3.^ (|3^^P at boundary AH (DE) is given. These angles are 
based on the input angles BETAI, /3j^, and BETAO, 

(4) These are calculated program constants, including the pitch from blade to blade, 
the mesh spacing, the minimum and maximum values of IT in the solution region 
(ITMIN and ITMAX), and the value of the prewhirl X (eq. (B8), ref. 2). 

(5) This is the number of mesh points in the entire solution region at which the 
stream function is unknown. 

(6) This is the boundary value (BV) of the stream function on each of the blade sur- 
faces. 

(7) This is the output corresponding to AANDK. 

(8) If the program calculates an optimum overrelaxation factor (i. e. , ORF = 0 in 
the input), then the successive estimates to the optimum value of ORF are printed. The 
last printed value of the estimated optimum ORF is the value of the overrelaxation factor 
(ORF) used by the program. 

(9) This is the output corresponding to ERSOR. 

(10) This is the output corresponding to STRFN. 

(11) This is the total execution time after obtaining the stream function solution for 
each outer iteration. 

(12) This is the output corresponding to SLCRD. 

(13) This is the output corresponding to INTVL. 

(14) This gives the maximum relative change in the density, for each outer iteration. 

(15) This is the output corresponding to SURVL. 
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(16) This is the total execution time after all calculations are completed for an outer 
iteration. 


ERROR CONDITIONS 

(1) SPLINT USED FOR EXTRAPOLATION 

EXTRAPOLATED VALUE = X.XXX 

SPLINT is normally used for interpolation, but may be used for extrapolation in some 
cases. When this occurs, the above message is printed as well as the input and output 
of SPLINT. Calculations proceed normally after this printout. 

(2) BLCD CALL NO. XX 

M COORDINATE IS NOT WITHIN BLADE 

This message is printed by subroutine BLCD if the m- coordinate given this subroutine as 
input is not within the bounds of the blade surface for which BLCD is called. The value 
of m and the blade- surface number are also printed when this happens. This may be 
caused by an error in the integer input items for the program. 

The location of the error in the main program is given by means of BLCD CALL 
NO. XX, which corresponds to locations noted by comment cards at each MHORIZ, 
ROOT, and BLCD call in the program. 

(3) ROOT CALL NO. XX 

ROOT HAS FAILED TO CONVERGE IN 1000 ITERATIONS 
This message is printed by subroutine ROOT if a root cannot be located. The input to 
ROOT is also printed. The user should thoroughly check the input to the main program. 

The location of the error in the main program is given by means of ROOT CALL 
NO. XX, which corresponds to locations noted by comment cards at each MHORIZ and 
ROOT call in the program. 

(4) DENSTY CALL NO. XX 

NER(l) = XX 

RHO*W IS X.XXXX TIMES THE MAXIMUM VALUE FOR RHO*W 
This message is printed if the value of pW at some mesh point is so large that there is 
no solution for the value of p and W. This indicates a locally supersonic condition, 
which can be eliminated by decreasing WTFL in the input. 

If RHO*W is too large, TURBLE still attempts to calculate a solution. This often 
permits an approximate solution to be obtained, which is valid at all the subsonic points 
in the region. In other cases the value of W is reduced at some of the points in question 
during later iterations, resulting in a valid final solution for these points. The program 
counts the number of times supersonic flow has been located at any point during a given 
run (NER(l)). When NER(l) = 50, the program is stopped. 
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The location of the error in the main program is given by means of DENSTY CALL 
NO. XX, which corresponds to locations noted by comment cards at each DENSTY call in 
the program, 

(5) MM, NBBI, NRSP, OR SOME SPLNO IS TOO LARGE 

If this is printed, reduce the appropriate inputs to their allotted maximum values. 

(6) WTFL IS TOO LARGE AT BLADE LEADING EDGE 

This is printed if WTFL is greater than the choking mass flow for the boundary BG. If 
this message is printed, WTFL is cut in half by the program and calculations proceed as 
usual. 

(7) ONE OF THE MH ARRAYS IS TOO LARGE 

This is printed if there are more than 100 intersections of horizontal mesh lines with any 
blade surface. In this case NBBI should be reduced. 

(8) THE NUMBER OF INTERIOR MESH POINTS EXCEEDS 2500 

This is printed if there are more than the allowable number of finite-difference grid 
points. Either MM or NBBI must be reduced. 

(9) SEARCH CANNOT FIND M IN THE MH ARRAY. 

If this is printed, the value of m and the blade-surface number are also printed. The 
user should thoroughly check the input to the main program. 


PROGRAM LISTING 


The program is identical to TANDEM (ref. 2) except for deleting all code dealing 
with the rear blade and making necessary corrections to the remainder. The program 
procedure for TANDEM, given in reference 2, is applicable also for TURBLE, except 
for small deletions. Also, the FORTRAN dictionary for TANDEM is valid for TURBLE 


COMMON SRW.ITER, I END , LER ( 2 ) , NER ( 2 J 

COMMON /AUKRHO/ A ( 2500, A ) , U ( 2500 » , K ( 2500) , RHO ( 2 500 ) 

COMMON /INP/GAM, AR, TIP, RHD IP, WTFL, OMEGA, ORE, BETA I.BET AO, 

1 MBI,MBO,MM,NB6I,NBL,NRSP,MR(50) ,RHSP(50) ,BESP(50>, 

2 BLDAT,AANDK,ERSOR,STRFN,SLCRD, INTVL.SURVL 

COMMON /CALC0N/M8IM1,MBIP1,MB0MI,MB0P1,MMM1,HM1 ,HT,DTLR,DMLR, 

1 PITCH,CP,EXPaN,TWW,CPTlP,TGROG,TBI,TBO,LAMBDA,TWL, ITMIN, ITMAX, 

2 NIP, IMS(2),BV(2),MV(100),IV(101),ITVI100,2),TV(100,2), 

3 DTOKVl 100,2) , BET AVI 100,2) ,MH 1 100, 2 ) , OTDMHI 100,2) ,BETAH( 100,2 ) , 

4 RMHI 100,2 ) ,BEH( 100,2 ) ,RMI 100) ,BE(100) ,DBDM( 100), SALIlOO), 

5 AAAI 100) 

COMMON /GEOMIN/ CHORD ( 2 ) ,STGR I 2 ) ,MLE 1 2 ) , THLE 1 2 ) ,RMH 2 ) ,RM0 ( 2 ) , 

1 RI12),R0I2),BETH2),BET0(2),NSPII2),MSPI50,2),THSPI50,2) 

COMMON /RHOS/RHOHBI 100,2),RHUVail00,2) 

COMMON /BLCDCM/ EM I 50, 2 ) , I NI T I 2 ) 

INTEGER BLOAT , A A NOK , ERSOR , STRFN , SLCRD, SURVL , A ATEMP, SORF , F I RST , 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,HV,MV1M1 
CALL TIMEliri) 
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10 lEND = -1 
ITER = 0 
iNirm = 0 
IN1T(2) * 0 
CALL MPUT 
CALL PRECAL 
30 CALL COEF 
CALL SDR 
CALL IIMEUT2) 

TIME= (T2-T1 )/3600. 

WRITEI6, 1000) TIME 
CALL SLAX 
CALL lANG 
CALL VELOCY 
CALL riMEl(T2) 

TIME= IT2-TD/3600. 

WRITE(6,1000) TIME 
IF(NERI2).GT.O) GO TO 10 
IF IIEMD) 30,30,10 

1000 FORMAT (8HLTIME = ,F7.4,5H MIM.) 
END 


SUBROUTINE INPUT 

input reads and prints all input data cards and CALCULATES HORIZONTAL 
SPACING IMV ARRAY) 

COMMON SRW, ITER, IEND,LER(2),NER(2) 

COMMON /AUKRHO/ A ( 2 500, A ) , U ( 2500 ) , K ( 2500 ) , RHOI 2500 ) 

COMMON /INP/GAM, AR,TIP,RHOIP,WTFL,OMEGA,ORF,BETAI,BETAO, 

1 MBI,MB0,MM,NBBI ,NBL,NRSP,MR(50),RMSPI50) ,BESP150) , 

2 bldat,aandk,ersor,strfn,slcro,intvl,survl 

COMMON /CALC0N/MBIM1,MSIP1 ,MB0M1 ,MB0P1 ,MMM1 ,HM1 ,HT, OTLR, DMLR , 

1 pitch,cp,expon,twm,cptip,tgrog,tbi,tbo,lambda,twl,itmin,itmax, 

2 NIP,IMS(2),BV(2),MV( 100), 1V( 101 ) , I T V ( 100, 2 ) , T VI 100, 2 ) , 

3 DTDMVI 100,2) , BET AV( 100,2) ,MH ( 100 , 2 ) , DTOMH 1 100 ,2 ) » BET AHI 100,2 ) , 

A rmHI 100, 2 ) ,BEH( 100,2 ) ,RM( 100 ) ,BE( 100) ,OBOM( 100) , SALIlOO), 

5 AAAI 100) , . 

COMMON /GEOMIN/ CHORD ( 2 ) , S TGR ( 2 ) ,MLE ( 2 ) , THLE ( 2 ) , RMl 1 2 ) , RMOI 2 ) , 

1 RI I 2) ,R0( 2 ),BETI ( 2), BETO( 2) ,NSPI (2) ,MSPI50,2 ) ,THSP(50,2) 

COMMON /RHOS/RHDHB (100,2) ,RHOVB ( 100,2 ) 

INTEGER BLDAT,AANDK,ERS0R,STRFN,SLCR0,SURVL,AATEMP, SURF, FIRST, 

1 UPPER, SI, ST, SRW 

REAL K,XAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
READ AND PRINT ALL INPUT DATA 
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WRITE(6,1000) 

READ( !>, UOO) 

WRITE(6, 1100) 

WRITE(6,1110) 

READ (5fl030) G AM , AR , T I P t RHO I P , WIFE » 8L ANK , OMEGA , ORF 
WRI TEI 6, 10^0) gam, AR,T I P,RHO IP, WTFL t BLANK, OMEGA, ORF 
WRITE(6, 1120) 

READ (5,1030)BETAI ,BETAO,CHORD( 1) ,STGR( 1) 

WR 1 1E( 6, 10^0) BET A I , BET AO, CHORD! 1 ) ,STGR! I ) 

WR1TE(6,1130) 

READ (5,1010) MBl ,MBQ, BLANK, BLANK, MM, N8BI,N8L,NRSP 
WRirE(6,1010) M8I ,M.80, BLANK, BLANK, MM, N08 1, NBL ,NRSP 
DO 10 J=l,2 

IF (J.EQ.l) WRITE!6, lUO) 

IF (J.eO.2) WRITE! 6, 1150) 

WRITE!6, 1180) J,J,J,J,J 

READ (5,1030) R I ! J ) , KO! J ) , BE T I ! J ) , BE TO ! J ) , SPLNO 
WRITE!6,1040) RI(J),RO!J),BFTI!J),BETO!J), GPLNO 
NSPI(J)= SPLNO 
NSP = NSPI ( J ) 

WR ITE ( 6, 1190) J 

READ (5,1030) ! M SP ( I , J ) , I = 1 , NSP ) 

WRITE!6, lOAO) !MSP( I , J) , I=1,NSP) 

WRITE(6,1200) J 

READ (5,1030) ( T HSP ! I , J ) , I = I , NSP ) 

10 WRITE(6, lOAO) ( THSP ( I , J ) , I = I ,NSP ) 

WR ITE(t>, 1210 ) 

READ (5,1030) ! MR ( I ) , 1= 1 ,NRSP) 

WRITE!5,1040) (MR( I ) ,I=1,NRSP) 

WRITE(L, 1220) 

READ (5,1030) ( RMSP ( I ) , I =1 ,NRSP ) 

WR ITE( 6, lOAO ) UMSP( 1 ), I=1,NRSP) 

WRITE(6,1230) 

READ (5,1030) ( B ESP ! I ) , 1 = 1 ,NRSP ) 

WR ITE( 6, 1040 ) ( BESP! I ) , I =1,NRSP) 

WRIIE(6,1240) 

READ (5,1010) BLOAT, AANDK,ERSOR,STRFN,SLCRD, INTVL.SURVL 
WR I IE! 6, 1020) BLOAT, AANOK,ERSOR,STRFN,SLCRO,INTVL,SURVL 
IF ( MM.LK. lOO.AND. NBB I .LE.50.ANO.NRSP.LE.50. AND.NSPI ( 1 ) .LE.50 
1 .AND. NSPI (2) .LE.50) GO TO 20 

WRITE (6,1250) 

STOP 

CALCULATE MV ARRAY 

20 HMl = CHORD! 1 )/FLOAT ( M3Q-MB1 ) 

DO 30 1M=1,MM 

30 MV(IM) = FLOAT! I M-MB I )*HM1 

C CALCULATE MISCELLANEOUS CONSTANTS 
C 

NER( 1)*0 
NER(2)=0 

PITCH = 2.*3.1415927/FLQAT(NBL) 

HT= PirCH/FLOAT! NBBI ) 

DTLR= HT/1000. 

DMLR = HMl/1000. 

BV!1) = 0. 

BV(2) = 1. 

M81Mi= MBI-1 
MB1P1= MBI+l 
MBOMl= MBQ-1 
MB0P1= M80+1 
MMMl = MM-1 
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CP = AK/ (GAM-1. ) *GAM 
EXPOM= l./(GAM-l.) 

TWW= 2.*0MEGA/WTFL 

CPTIP= 2.*CP*TtP 

TGRUG= 2.*GAM*Afl/( GAM-K . ) 

CALL SPL IMT(MR,RMSPtNRSP|MV,MM,RMtSAL) 

CALL SPUNK MR.BESPtNRSPfMV.MM.BEtOBDM) 

CALCULATE GEOMETRICAL CONSTANTS 

CHORD(2 I = CHORD(l) 

STGR12) = STGR( 1 ) 

MLE( 11=0. 

MLE(2) = 0. 

THLL(l) = 0. 

THLEC2) = PITCH 
RMIll) = RM(MBt) 

RMK2) = RM(MBI) 

RMU( 1) = RM(MBO) 

RM0(2) = RM(MBO) 

INITIALIZE ARRAYS 

DO 60 1=1,2500 
0(1) = I. 

K( I 1 = 0. 

60 RHO( I ) = RHOIP 
00 70 IM=1,100 
DO 70 SURF=1,2 
RHOHBI IM.SURF ) = RHOIP 
70 RHOVBI IM.SURF ) = RHOIP 
RETURN 

1000 FORMAT (IHl) 

1010 FORMAT 11615) 

1020 FORMAT ( IX, 1617 ) 

1030 FORMAT (8F10.5) 
lOAO FORMAT (1X,8G16.7) 

1100 FORMAT (80H 

1 ) 

1110 FORMAT ( 7X , 3HGAM ,1AX,2HAR, 13X,3HTIP, 12X,5HRH01P, 12X , AHWTFL , 1 IX , 6HW 
1TFLSP,10X,5H0MEGA, 12X,3H0RF ) 

1120 FORMAT ( 6X , 5H8E T A I , 1 OX, 5HBET AO , 1 IX , 6HCH0R0F , 1 1 X , 5HSTGRF ) 

1130 FORMAT (<,10 MB I MBO MM NBBI NBL NRSP) 

llAO FORMAT { 39HL BLADE SURFACE I — UPPER SURFACE) 

1150 FORMAT { 39HL BLADE SURFACE 2 — LOWER SURFACE) 

1180 FORMAT (7X,2HRI, I1,12X,2HR0,I1, 12X,',HBET I , 1 1 , 1 1 X ,AHBETO, I1,11X,5HS 
IPLNO, II > 

1190 FORMAT ( 7X , 3HMSP , I I , 2X, 5HARRAY ) 

1200 FORMAT ( 7X , 4HTHS P , I 1 , 2X , 5HARR AY ) 

1210 FORMAT ( 16HL MR ARRAY) 

1220 FORMAT (7X,11HRMSP ARRAY) 

1230 FORMAT (7X,11HBESP ARRAY) 

I2A0 FORMAT ( 52HL BLDAT AANDK ERSOR STRFN SLCRD INTVL SURVL ) 
1250 FORMAT lAlHl MM , NBH 1 , N 1 SP , OR SOME SPLNO IS TOO LARGE) 

END 
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SUBROUriNE PRECAL 
C 

C PRECAL CALCULATES ALL REQUIRED FIXED CONSTANTS 
C 

COMMON SRW, ITER, I END , LER ( 2 » , NER I 2 I 

COMMON /AUKRHO/ A ( 2 500, 4 ) , U< 2500 ) , K ( 2500 1 , RHO ( 2 500 ) 

COMMON /INP/GAM, AR,TIP,RH0IP,WTFL,0MFGA,0RF,BETAI,3ETAU, 

1 MBI , MBO,MM,NBBI , NBL, NRSP,MR( 50) ,RMSP( 50) ,6ESP(50) , 

2 BLOAT , AANDK,ERSQR,STRFN,SLCRD, INTVL.SURVL 

COMMON /CALC0N/MBIM1,MBIP1,MB0MI,MB0P1 ,MMMI,HM1 ,HT , DTLR , OMLR , 

1 P1TCH,CP,EXPON,TWW,CPTIP,TGROG,TBI,TBO,LAMBDA,TWL, ITMIN, ITMAX, 

2 NIP, IMS( 2 ) ,BV (2) ,HV( 100), IV( 101 ) , ITV( 100,2) ,TV( 100,2) , 

3 DTDMVI 100,2), BET A VI 100,2 ) ,MH« 100,2) ,DTDMH( 100,2 ) ,BET AH ( 100,2 ) , 

4 RMH( 100,2), BEH( 100,2) ,RMt 100) ,BE( 100) ,0B0M(100) , SAUlOO), 

5 AAAI 100) 

INTEGER BLOAT , AANDK , ERSOR , STRF N , SLCRD, SURVL , AAT EMP, SURF , F I RST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLF,MR,MSL,HSP,MV,HVIM1 
EXTERNAL BL1,BL2 
C 

C CALCULATE LAMBDA AND VI 
C 

BETAI = BETAI/57. 295779 
BETAO = 8ETA0/57. 295779 
TBI =« S INI BETAI ) /COS (BETAI ) 

TBO = SINI BETAO) /COS I8ETA0) 

10 RHOT = RHOIP 

RHOVI = WTFL/BEIMBI )/PITCH/COS(BETAI )/RH(MBI ) 

20 VI = RHOVI/RHOT 

LAMBDA = RM(MBI)*IVI*SINIBETAI)+0MEGA*RM(MB1)) 

TTIP = 1.- ( V I**2+2.*OMEGA*LAMBOA-IOMEGA*RM (MBI ) )**2 ) /CPTIP 
IFITTIP.LE.O. ) GO TO 30 
RHOMBI = RHOIP*TTIP**EXPON 

IF(ABS(RH0MBI-RH0T)/RH3IP.LT. .000001) GO TO 40 
RHOT = RHOMBI 
GO TO 20 

30 MTFL = WTFL/2. 

NER(2)= NERI2)+l 
MR1T£I6,1020) WTFL 
IFIN£RI2).EQ.10) STOP 
GO TO 10 

40 VI = RHOVI/RHOMBI 

LAMBOA = RMIMBI ) ♦( VI4S1NIBETAI )+OMEGA*RM(HBl ) ) 

C 

C CALCULATE MAXIMUM VALUES FOR RHO*H AT LEADING AND TRAILING EDGE 
C 

TWL = 2.*0MEGA*LAM3DA 

AA = (TWL-(OMEGA*RM(M0I ) )**2)/CPTIP 

TPP = TIP*! l.-AA ) 

BB = TGROG*TPP 
TTIP = l.-BB/CPTIP-AA 
RHOT =1 RHOIP*TTI P**EXP3N 
RHOWMI = RHOT*SQRT(BB ) 

AA = (TWL-(0MEGA*RM(MB3 ))**2)/CPTIP 
TPP = TIP*I1.-AA) 

BB = TGROG*TPP 
TTIP = l.-BB/CPTIP-AA 
RHOT = RH0IP*TTIP**EXP3N 
RHOWMO = RHOT*SORT( BB) 

CALCULATE VO AND W-CRITICAL AT BLADE LEADING AND TRAILING EDGE 
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RHOVO = WTFL/BE(MBO ) /P I TCH/COS < BE TAO) /RM( MBO » 

RH0HB2 = RHOIP 

TWLMR = TWL-(OMEGA*RM(MBO ))**2 
LER( l)=l 

C DENSTY CALL NO. 1 

CALL OENSTYI RHOVO, RHOMB 2 , VO, TWLMR ,CPTI P,EXPON, RHOI P , GAM, AR, T 1 P ) 
WCRI = SQRT( TGR0G*TIP*( l.-(TWL-(OMEGA*RM(M8I ) »**2)/CPnP) ) 

WCRO ^ SQRT( TGROG*TIP*( 1.- (TWL-(OMEGA*RM(MBO ) ) **2 ) /CPT I P ) ) 

CALCULATE BETA CORRECTED TO BOUNDARY A-N AND G-H 

TWLMR =« TWL-(OMEGA*RM(l ) )*«2 
RHOl = RHOMB! 

TBIl = 1.E20 

50 TBIT = (TBI/BE(MBI )*RHDl/RH0MBI»-0MEGA*(RM(MBn**2-RM( I (♦*2)*RHU1 
1 /WTFL*P1TCH)*BE( I ) 

IF(ABS(TBI1-TBIT).LT. .00001) GO TO 60 
TBIl = TBIT 

RHOVI = WTFL/PITCH*SQRT( l.+TBIl**2)/BE(l)/RM(l) 

LER( U=2 

DENSTY CALL NO. 2 

CALL DENSTY ( RHO VI , RHO I , AA , T WLMR ,C PT I P ,EXPON, RHOI P, GAM , AR , T 1 P ) 

GO TO 50 
60 TBl = TBIT 

STAIN = ATAN( TBI )*57. 295779 
TWLMR = TWL-(OMEGA*RM(MM) )**2 
RHOMM = RH0MB2 
TBOM = I.E20 

70 TBOT = (TB0/8E(M80 )*RHOMM/RHOMB2+QMEGA*(RM(MBO I ♦*2-RMI MM) ) ♦ 
I KHOMM/WTFL*PI TCH)*BE(MM) 

IF (ABS(TBOM-TBOT).LT. .00001) GO TO 80 
TBOM = TBOT 

RHOVO = WTFL/PITCH*SORT( l.+TB0M»*2)/8E(MM)/RM(MM) 

LER! 1)=3 

DENSTY CALL NO. 3 

CALL DENSTY ( RHO VO, RHOMM , A A, TWLMR , CPT I P, EXPON , RHOI P .GAM, AR.TIP) 
GO TO TO 
80 TBO = TBOT 

BTAOUT = ATAN(TB0)*57. 295779 

CALCULATE TV, ITV, IV, DTDMV, AND BETAV ARRAYS 

. ITMIN = 0 

ITMAX = NBBI-1 
C TV, ITV, AND DTDMV ON BLADE 
DO 90 IM=MBI,MBO 
LERI2)=l 

C BLCD CALL NO. 1 

CALL BLUMVI IM),TV( IM,1 ) , DTDMV! IM,l),INF) 

ITVI IM, 1 )= INTI ( TV( IM,l )*DTLR)/HT) 

IF ITVI IM, D.GT.-OTLR) I TV 1 1 M, I )» I TV ( IM, I ) ♦ I 
ITMIN= MINOI ITMIN, ITVI IM, I ) ) 

LER(2)=2 

C BLCD CALL NO. 2 

CALL BL2(MV(IM), T V II M , 2 ) , DTDMV ( I M, 2 ) , I NF ) 

ITVI IM,2)= INTI I TV! I M, 2 ) -OTLR ) /HT ) 

IF ITVI IM,2).LT.DTLR) I T V I I M , 2 ) I T VI I M, 2 ) - 1 
90 ITMAX= MAX0I1TMAX,ITVIIM,2)) 

C ITV AND IV UPSTREAM OF BLADE 
FIRST = 0 
LAST = NBBI-l 
DO 120 IM*1,MB1M1 
ITVI IM, 1 )= FIRST 
120 ITV(IM,2)= LAST 
C ITV DOWNSTREAM OF BLADE 
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140 LAST= ITVJMB0,2) 

FIKST= LAST+l-NBBI 
DO 150 IM=M80Pl,MM 
1TV( IM, I ) = FIRST 
150 ITV( IM.2) = LAST 

ITMIN = MINOI ITMIN, ITV( MM, in 
C CALCULATE IV ARRAY 
1V( 1 ) = 1 
00 160 IM=1,MM 

160 IVIIM + ll = IV(IM) + ITV(lM,21-ITV(IM,im 
C BETAV ARRAY 

00 200 SURF=1,2 
DO 200 IM=M8I,M80 

200 BETAV( IM.SURF) = AT AN ( OT DMV I IM, SURF I ♦RMI I M > ) *57 . 295779 
NIP = IV(MM)+NBBI-1 

WRITE! 6, 10 30) VI .RHOWMI , WCRI , 8T A I N , VO, RHOWKO, MCRO, BT AOUT 

WRITEI6, 1040) PITCH, HT.HMl 

WRITE (6, 1050) ITMIN, 1 TM AX , LAMBDA , N IP 

WRITE (6, 1060) ( SURF, RV( SURF) ,SURF= 1,2) 

IF(BLOAT.LE.O) GO TD 230 
WRITE (6,1070) 

WRITE (6,1080) (MV( IM),TV( IM,1) ,DTDMV( I M , 1 ) , TV ( I M, 2 ) , DT DMV (I M, 2 ) , 
1 IM=MBI,M80) 

WRITE (6,1090) (IM,MV(IM),RM(IH), SAL( IM) ,BE (I M ) , DBDM ( I M ) , I M= 1 , MM ) 
230 CONTINUE 
C 

C CALCULATE HH AND DTDMH ARRAYS 
C 

ITO = I TV( I, 1 ) 

MRTS = I 
IMS(l) = 1 
MH(l,l) = 0. 

DTDMH(1,1) = I.EIO 
LER(2) = 3 

C BLCD AND ROOT (VIA MHORIZ) CALL NO. 3 

CALL MH0RIZ(MV,ITV(1,1),BL1,M6I,MB0,IT0,HT,0TLR,0,IMS(1),MH(1,1), 
1 DTDMHI 1 , 1 ) ,MRTS ) 

IF ( ITV(MBO,l)-ITV(MBO,2)+NBBI.NE.2) GO TO 240 
IMSL = IMS( I)+l 
MH(IMSL,1) = MV(MBO) 

DTDMHI IMSL, 1 ) = -l.ElO 
IMS(l) = IMSL 
240 IMS(2) = 0 
MRTS = 1 
LER(2) = 4 

C BLCD AND ROOT (VIA MHORIZ) CALL NO. 4 

CALL MHORIZ (MV, I TV ( 1,2) , BL2 , MB I , MBO, I TO, HT ,OTLR , 1 , I MS I 2 ) , MH ( 1, 2 ) , 

1 DTDMHI 1 ,2 ) ,MRTS ) 

I = MAX0( IMS( 1), IMSI 2) ) 

IF( I.LE.lOO) GO TO 290 

WRITE(6,1100) I 

STOP 

290 IF(BLDAT.LE.O) GO TO 300 

WRITE (6,1110) ( IM, I V( IM) , ( ITVI IM,SURF ) ,SURF=1 ,2) , 1M= 1 , MM ) 

CALCULATE RMH, BEH, AND BETAH ARRAYS 

300 IF(BLOAT.GT.O) WR I TE ( 6, 1 120) 

DO 320 SURF=1,2 

CALL SPLINT! MR, RMSP, NR SPfMHIl, SURF), IMSI SURF ) ,RMH(l, SURF), AAA) 
CALL SPLINT(MR,BESP,NRSP,MH( 1 , SURF ), IMSI SURF),BEH(1, SURF), AAA) 
IMSS = IMSISURF) 

IF I IMSS.LT. 1 ) GO TO 320 
DO 310 IHS = 1,IMSS 
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310 8ETAH( IHSfSURF) = AT AN I DTOMH ( IHS . SURF ) *RMH ( I HSt SURF ) ) *5 7 .295 779 
IF (BLOAT. GT.O) HR I TE ( 6 , 1 1 30 ) SURF , ( MH t I M , SURF ) , RMH (I M , SURF ) , 

1 BEH( IMf SURF) , BE T AH( I M t SURF ) , OTOMH (I M , SURF ) , IM^lt IMSS) 

320 CONTINUE 

IF (BLOAT. LE.O) GO TO 340 
WRITE (6,1140) 

IT = ITMIN 

330 IF ( IT.GT. ITMAX) GO TO 340 
TH = FLOAT ( IT )*HT 
WRITE (6,1010) IT,TH 
IT = IT+1 
GO TO 330 

340 IF (NIP. LE. 2500) GO TO 350 
WRITE (6, 1150) 

STOP 

350 WRITE (6,1000) 

RETURN 

1000 FORMAT (IHl) 

1010 FORMAT (4X, I4,G16.5) 

1020 FORMAT! 60HL INPUT WEIGHT FLOW (WTFL) IS TOO LARGE AT BLADE LEADING 
1EDGE/16H WTFL REDUCED T0,G14.6) 

1030 FORMAT ( 1H1/24X , lOHFREESTREAH, 8X, 1 3HHAXIMUM VALUE, 

17X,faHCR ITICAL, 30X, 14HBETA CORREC TED/25X , 8H VELOC I TY, lOX , 9HF0R RHO*H 
2, lOX.BHVELOCITY, 31X, IIHTO BOUNOARY/1 X, 17HLEADING EDGE B-G,3G18.5, 
312X, 12HB0UNDARY A-H, G 18 . 5/ IX , 1 7HTRAI LI NG EDGE C-F, 3G 18 . 5 , 12X , 
412HBUUNDARY D-F,G18.5) 

1040 FORMAT! 33HL CALCULATED PROGRAM CONS TANTS//5X , 5HP I TC H , 1 3X , 

1 2HHT,13X,3HHH1/1X,5G16.7I 

1050 FORMAT (/5X,5H ITMIN, lOX , 5H I TMAX/4X , I 5 , lOX, 1 5//5 X ,6HL AMRDA/ 1 X , G 1 6. 7 
1 ,/38HL NUMBER OF INTERIOR MESH POINTS = ,15) 

1060 FORMAT! 28HL SURFACE BOUNDARY VALUES// 5X , 7HSURF ACE , 7X , 2HBV 
1/(5X, I4,4X,F10.5 ) ) 

1070 FORMAT ( IH 1 , 6X, 6 2HBL ADE DATA AT INTERSECTIONS OF VERTICAL MESH LIN 
lES WITH BLADES) 

1080 FORMAT ( IHL , 22X , 15HBLA0E SURFACE 1 , 15X , 15HBL ADE SURFACE 2/7X, 

1 IHM, 14X.2HTV, 11X,5HDT0MV, 12X , 2HT V , 1 IX , 5HDTDMV/ ( 5G 15 . 5 ) ) 

1090 FORMAT ( IH 1 , 1 3X , 44HS TRE AM SHEET COORDINATES AND THICKNESS TABLE / 

1 2X,2HIM,7X,lHM,14X,lHR,13X,3HSAL,13X,lHB,12X,5HDB/DM/( IX, 13, 

2 5G15.5)) 

1100 FORMAT ( 34HL0NE OF THE MH ARRAYS IS TOO LARGE/7H IT HAS, 15, 8H POI 
INTS) 

1110 FORMAT (4H1 IM,9X,8HIV ARR AY ,2 5X , 9HI TV ARRAY/ 38X , 5HBL ADE/ 37X , 7HSUR 
IFACE, 3X, 1H1,5X, 1H2/39X, 3HN0. / U X, I 3, 5X , 1 10 ,25X , 2 ( 14 , 2X ) ) ) 

1120 FORMAT (67H1M COORDINATES OF INTERSECTIONS OF HORIZONTAL MESH LINE 
IS WITH BLADE) 

1130 FORMAT (25HLMH ARRAY - BLADE SURFACE , 1 2/ / 1 5X , 2HMH, 19X, 3HRMH , 19X , 

1 3HBEH,18X,5HBETAH,17X,5HDTDMH/(5G22.4) ) 

1140 FORMAT (43H1THETA COORDINATES OF HORIZONTAL MESH L I NES// 6X , 2H I T , 
15X,5HTHETA) 

1150 FORMAT! 48HLTHE NUMBER OF INTERIOR MESH POINTS EXCEEDS 2500) 

END 
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SUBROUTINE MHOR I Z ( NV , I T V tBL , MB I ,MBO, I TO , HT ,OTLR ,KODE , J , MH, DTDMH, 
IMRTS ) 

C 

C MHORIZ CALCULATES M COORDINATES OF INTERSECTIONS OF ALL HORIZONTAL 
C MESH LINES WITH A BLADE SURFACE 
C KUDE = 0 TOR UPPER BLADE SURFACE 
C KUDE =* I FOR LOWER BLADE SURFACE 
C 

COMMON SRW.ITER, I END, LER ( 2 I , NER ( 2 I 
DIMENSION MV( 100 ) , I TVl TOO) ,MH( 100) ,OTDMHIlOO) 

INTEGER BLOAT,AANDK,ERSOR,STRFN,SLCRD, SURVL, AATEMP, SURF, FIRST, 

I UPPFR,S1,ST,SKW 

REAL K,KAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,MV,MV1M1 
REAL MVIM 
EXTERNAL BL 

IF (MBI.GE.MBO) RETURN 
IM= M81 
10 ITINO= 0 

20 IF (ITVl IM+l )-ITV( IM)-I TIND) 30, <10, 50 
30 J= J+1 

TI= FLOAT! 1TV(IM+1)-ITQ-ITIN0«^K0DE)*HT 
ITIND= ITINO-1 
MVIM = MV( IM ) 

IF (MKTS-EQ.1) MVIM = MVIM+IMVI IM+l)-MVIM)/lOOO. 

CALL ROOT (MVIM,MV(IM^l),TI,BL,DTLR,MH(J),DTDMH(jn 
60 TO 20 
AO IM= IM+l 
MRTS = 0 

IF (IM.EQ.MBO) RETURN 
GO TO 10 
50 J= J+1 

TI= FLOAT! 1TVIIM)-IT0+ITIND+K00E)*HT 
ITIN0= ITINO+1 
MVIM = MV! IM) 

IF IMRTS. EQ.l) MVIM = MVIM+! MV! IH+1 )-MVIM) /lOOO. 

CALL ROOTIMVIM , MV I IM+l ) , T I , BL , OTLR ,MH ! J ) , DTOMH ! J ) ) 

GO TO 20 
END 
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$I8FTC COEF 


DEBUG 


SUBROUTINE COEF 

COEF CALCULATES FINITE DIFFERENCE COEFFICIENTS, A, AND CONSTANTS, K, 
AT ALL UNKNOWN MESH POINTS FOR THE ENTIRE REGION 

COMMON SRW, ITER, IEND,LER(2),NER(2) 

COMMON /AUKRHO/ A ( 2500, A ) , U( 2500 ) , K ( 2500 ) , RHOI 2500 I 
COMMON /INP/GAM, AR , T I P, RHOI P , WTFL, OMEGA , ORF , BE T A I , BETAC, 

1 MBI ,MBO,MM,NBBI ,NBL,NRSP,MR(50) ,RMSP( 50) ,BESP(50) , 

2 BLUAT,AANOK,ERSOR,STRFN,SLCRD,INTVL,SURVL 

COMMON /CALC0N/MBIM1,MBIP1,MB0M1,MB0PI,MMM1,HM1 ,HT,DTLR,DMLR, 

1 PITCH,CP,EXPON,TWW,CPTIP,TGROG,TBI,TBO,LAMBDA,TWL,ITMIN, UMAX, 

2 NIP, IMS(2) ,BV 12) ,MV( 100) , IV( lOl) , IT V ( 100,2 ) , T VCl 00 , 2 ) , 

3 DTOMVI 100,2), BETAVI 100,2) ,MH( 100,2) , DTDMH ( 100 ,2 ) , BE T AH( 100 , 2 ) , 
A RMH( 100,2) ,BEH( 100,2 ) ,RM( 100) ,BE( 100) ,OBDMl 100) , SAL(IOO), 

5 AAAI 100) 

COMMON /HRBAAK/H (A) ,RIA) ,B(A),KAK(A) ,KA I A) ,RZ , B Z , IH ( A ) 

INTEGER BLOAT, AANDX,ERSOR,STRF.N,SLCRD,SURVL,AATEMP, SURF, FIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMB0A,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
INITIALIZE ARRAYS 
ITER = ITER+l 
IHIl) = I 
IHI2) = 0 

INCOMPRESSIBLE CASE 

IFIGAM.NE. 1.5.0R.AR.NE. 1000..0R.TIP.NE.I.E6) GO TO 20 
lEND = 1 
GO TO AO 

ADJUSTMENT OF PRINTING CONTROL VARIABLES 
20 1F(1TER.NE.1.AN0.ITER.NE.2) GO TO 30 
AANDK = AANDK-l 
ERSOR = ERSOR-1 
STRFN = STRFN-l 
SLCRD = SLCRD-l 
INTVL = INTVL-1 
SURVL = SURVL-1 
30 1F( lENO.NE.O ) GO TO AO 
AANDK = AANDK+2 
ERSOR = ERSQR+2 
STRFN = STRFN+2 
SLCRD = SLCRD+2 
INTVL = INTVL+2 
SURVL = SURVL+2 

FIRST VERTICAL MESH LINE 

AO DO 50 IP=1,NBBI 
A( IP,1) = 0. 

A( IP, 2) = 0. 

AI IP, 3) = 0. 
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AdPf'*! = 1. 

50 K(IPJ = HM1*TB1/PITCH/RM( 1 ) 

UPSTREAM OF BLADE, EXCEPT FDR FIRST VERTICAL MESH LINE 

IF (2-GT.MBIMl ) GO TO 70 
DO 60 IM=2,MBIMi 
60 CALL C3EFP( IM) 

BETWEEN BLADES 

70 DO 80 IM=MBI,MB0 
80 CALL CDEF88( IM) 

DOWNSTREAM OF BLADES EXCEPT FOR FINAL MESH LINE 

150 IFIMBOPl.GT.MMMl ) GO TO 170 
DO 160 IM=MB0P1, MMMl 
160 CALL COEFPl IM ) 

FINAL VERTICAL MESH LINE 

170 IVMM = IV(MM) 

DO 180 IP=IVMM,NIP 
A( IP,1) = 0. 

A( IP, 2) = 0. 

A( IP, 3) = 1. 

A( IP, A) = 0. 

180 K(IP) = -HMl*TBO/PITCH/RM(MM) 

TAKE CARE OF POINTS ADJACENT TO B, AND CASES WHEN POINTS J,C,E, OR F 
ARE GRID POINTS 

POINT B 

IP = IVIMBIMI ) 

A(IP,4) = 0. 

C POINT C 

IF( ITV(MBO,1)-ITV(MBO,2)+NBBI.NE.2) RETURN 
IT = ITVIMBO, l)-l 
IP = IPFIMBOPl, I T I 
A ( I P , 3 ) = 0 . 

RETURN 

END 
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tIBFTC COEFBB DEBUG 

SUBKOUriNE CDEFBB( IM) 

COEFBB CALCULATES FINITE DIFFEKENCE COEFF I C lENTS t A, AND CONSTANTS, 
ALONG ALL VERTICAL MESH LINES WHICH INTERSECT BLADES 

COMMON /AUKRHD/ A ( 2500, A ) , U ( 2500 ) , K I 2500 ) , RHO I 2 500 ) 

COMMON /INP/GAM, AR, TIP, RHO IP, WTFL, OMEGA, ORF, BETA 1, BE TAG, 

1 MBI,MB0,MM,NB8I,NBL,NRSP,MR( 50) ,RMSPI 50 ) , BE SP I 50 ) , 

2 BLDAT,AANOK,ERSOR,STRFN,SLCRO,INTVL,SURVL 
COMMON /CALCON/MBIMl.MB I P 1 ,MBOM I , MBQPl , MMMl , HMl , HT , OTL R , DMLR , 

1 PITCH,CP,EXPON,TWW,CPTIP, TGROG, TB I , TBO, LAMBDA , TWL , I IMIN, ITMAX, 

2 NIP, IMS( 2 ) ,BV(2) ,MV( 100) , IVI 101 ) , ITVI100,2) ,TV( 100,2) , 

3 DTDMVI 100, 2) , BETA V( 100,2) ,MH 1100, 2) ,OTDMHI 100 ,2 ) , BE T AH ( 100 , 2 ) , 
A RMH( 100,2) ,BEH( 100,2 ) ,RM( 100) ,BE( 100) ,OBDM( 100) , SAL(IOO), 

5 AAAI 100 ) 

COMMUN /HRBAAK/HI4) ,Rt A) ,B(A),KAK(A) ,KA(A),RZ,BZ,IH(A) 

INTEGER 8LDAT,AANDK,ERS0R, STRFN,SLCRO,SURVL,AATEMP, SURF, FIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMB0A,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
IF( 1 TV! IM, 1 ) ,GT. ITVI IM, 2 ) ) RETURN 
ITVU = ITVI IM,1) 

IT = ITVU - 1 
ITVL = ITVI IM,2) 

IPU = IPFI IM, ITVU) 

IPL = IPU+ITVL-ITVU 
DO 90 1P=IPU,IPL 
IT = IT+1 

CALL HRBI IM, IT, IP) 

DO 10 1=1, A 
KAKIl) = 0. 

10 KAI I) = 0 

C FIX HRB VALUES FOR CASES WHERE MESH LINES INTERSECT BLADES 
60 IFIIT.EO.ITVIIM, in CALL BOR Y12 ( 1 , I M, I T ) 

IFIIT.EQ.ITVI IM,2) ) CALL BOR Y12 I 2 , I M , 1 T ) 

ITVMl = ITVI IM-1 , 1 ) 

ITVPl = ITVIIM+1,1) 

IFI1T.lt. ITVM l) CALL BDR Y3AI 3, I M, 1 ) 

IFllT.LT. ITVPl) CALL B0RY3AI A, I M, 1 ) 

IFI IT.GT.ITVI IM-1,2) ) CALL BDRY3A1 3, IM,2 ) 

IFI IT.GT.ITVI IM<-1,2) ) CALL BDRY3AIA, IM,2 ) 

70 IFI IM.EQ.MB0.AND.LQWER.EQ.2) GO TO 80 
C COMPUTE A AND K CDEFFICIENTS 
80 CALL AAKI IM, TP) 

DO 90 1=1, A 

KIIP) = KI IP ) +KAKI I ) *AI IP, I ) 

90 IFIKAI I ) .EQ. 1 ) AIIP,I) = 0. 

RETURN 

C 

C COEFP CALCULATES FINITE DIFFERENCE COEFFICIENTS, A, AND CONSTANTS, K 
C ALONG ALL VERTICAL MESH LINES WHICH DO NOT INTERSECT BLADES 
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ENTPY CGEFP( IM) 

ITVU = ITVdM.l) 

IT = lTVU-1 
ITVL = 1TV(IM,2) 

IPL = I V( IM + 1 )-l 
IPU = I VUM) 

DO 100 lP=IPUfIPL 
IT = II+l 

CALL H:T8( IM, IT, I PI 
IF (IT. EQ. ITVU) Ptl) = RHO(IPL) 
IF (IT. EQ. ITVL) R(2) = RHO(IPU) 
100 CALL AAK( IM, IP) 

K( IPL ) = K( IPL) + A( IPL,2 ) 

K(IPU) = K( IPU)-A( IPU,1 ) 

RETURN 

END 


tIBFTC HRH DEBUG 

SUBROUTINE HRB( IM, IT,!P) 

HRB CALCULATES MESH SPACING, H, DENSITIES, R2 AND R, AT GIVEN AND 
ADJACENT POINTS, AND STREAM SHEET THICKNESSES, 8Z AND B, AT GIVEN 
AND ADJACENT POINTS 

COMMON /AUKRHO/ A ( 2 500, A ) , U( 2500 ) , K I 2500 ) , RHO ( 2500 ) 

COMMON /CALC0N/MBIM1,MBIPI ,MBQMl ,M80P1,MMMI,HM1 ,HT,DTLR,DMLR, 

1 PITCH,CP,EXPON,TWW,CPTlP, TGROG, TB I , TBO, L AMBDA ,TWL , 1 TMI iN, I TM AX 

2 NIP,IMS(2),8V(2),MV(100),IV(101),ITV( 100,2 ) , TV( 100,2 ) , 

3 OTDMVI 100, 2) , BET AVI 100,2) ,MH( 100,2) ,DTDMHI 100 ,2 ) , BET AH ( 100 , 2 ) 
A KMH( 100,2 ) ,BEH( 100,2 ) ,RM( 100), BE(IOO) ,0BDMI100) , SALIlOO), 

5 AAAI 100) 

COMMUN /HRBAAK/H ( A),R(A),B(AI,KAK(A),KA(AI,RZ,BZ,IH(A) 

INTEGER BLDAT , A ANDK , ERSOR , STRF N , SLCRO , SURVL , AAT EMP, SURF , F I RS T , 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
H( 1)= HT*RM( IM) 

H( 2)= HT*RM( IM) 

HI 3)= MVI IM) - MVl IM-l) 

H(A) = MVI IM+D-MVI IM) 

RZ = RHOIIP) 

IP3 = IPFI IM-l, IT) 

IPA = IPFI IM+1, I T ) 

RID = RHQI IP-1) 

RI2) = RHOI IP*1) 

Rl 3) = RHOI IP3) 

Rl A) = RHOI IPA) 

BZ= 8EI IM) 

Bl 3) = HE I IM-l ) 

BIA)= BEI IM + 1 ) 

RETURN 
END 



$IBFTC AAK 


DEBUG 


SUBKOUrtME AAK( IM, IP) 

C 

C AAK CALCULATES FINITE DIFFERENCE COEFFICIENTS, A, AND CONSTANT, K, 

C AT A SINGLE MESH POINT 
C 

COMMON /AUKRHO/ A ( 2500 , A ) , U ( 2500 ) , K ( 2500 ) , RHO 1 2 500 ) 

COMMON /CALCON/MBIMI,MBIP1,MBOMI,MBOP1 ,MMMl ,HM1 ,HT ,OTLR,DMLR , 

1 P1TCH,CP,EXP0N,TWW,CPT1 P,TGROG,TBI , TBG, L AMBDA ,TWL , I TM I N, I TMAX , 

2 NIP, IMS! 2 ) ,B\M2 ) ,MV( 100) , IV ( 101 ), ITVI 100,2 ) ,TV( 100,2) , 

3 DTDMVI 100,2) , BET A V ( 1 00 , 2 ) , MH ( 100, 2 ) , 0 TOMH ( 100 ,2 ) , BET AH { 1 00 , 2 ) , 

4 RHHI 100, 2 ) ,BEH( 100,2 ) ,RM( 100) ,BEI 100) tOBDMI 100) , SAL(IOO), 

5 AAAI 100) 

COMMON /HRBAAK/HIA) ,R(',) , B 1 4 ) , KAK ( 4 ) , K A ( 4) ,RZ,BZ, IH(^) 

INTEGER BL DAT, A ANDK, ERSOR, STRFN,SLCRO,SUHVL, AATEMP, SURF, FIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMB0A,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
A12= 2. /HI D/H12) 

A3A= 2./HI 3)/HU) 

AZ= A12+A3^ 

B12= I Rl 2 )-RI 1) ) /RZ/ IHI 1 )+Hl 2) ) 

B3A= IBIA)*RIA)-BI3)*RI3)) /B Z/RZ/ I H I 3 ) +H I A ) ) -SAL I I M ) /RM I I M ) 
AIIP,1) = I2./HI l)4^ei2)/AZ/IHIl)+HI2)) 

AIIP,2) = A12/AZ-AI IP,1 ) 

AIIP,3) = I 2./HI 3)+B34)/AZ/IH( 3I+HIA) ) 

AI IP, A) = A3A/AZ-AI IP,3 ) 

KIIP) = -TWW*BZ*RZ*SALI IM)/AZ 

RETURN 

END 


$IBFTC BDRY12 DEBUG 

SUBROUTINE B DRY 1 2 I I , I M, I T ) 

C 

L BDRY12 CORRECTS VALUES COMPUTED BY HRB WHEN A VERTICAL MESH LINE 
C INTERSECTS A BLADE 
C 

COMMON /CALC0N/MBIM1,M3IP1 ,MB0M1 ,MB0P1 ,MMH1 , HMl ,HT , DTLR, OMLR , 

1 PITCH,CP, EXPON, TWW,CPTI P, TGROG , TB I , TBO, LAMBDA , TWL , ITMIN,ITMAX, 

2 NIP, IMSI2 ) ,8V 12) ,MVI 100 ) , I V I 101 ) , I TV I 100, 2 ) , T VI 1 00 , 2 ) , 

3 DTDMVI 100, 2), BET AVI 100,2 ) ,MHI 100,2) .OTDMHI 100,2 ) .BETAHI 100,2 ) , 
A RKHI 100,2),BEHI100,2),RMI 100),eEI100),D60MI100), SAL(IOO), 

5 AAAI 100) 

COMMON /RHOS/RHOHBI 100, 2 ) ,RHDV8 I 100,2) 

COMMON /HRBAAK/H I A) ,RIA),BIA),KAKIA),KAIA),RZ,BZ,THIA) 

INTEGER BL 0 AT, A A NDK, ERSOR, STRFN,SLC«0,SURVL, AATEMP, SURF, FIRST, 

I UPPER, SI, ST, SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
HID = ABSIFLOATI I T)*HT-TVI IM,I n*RMITM) 

RI I )= RHOVBI IM, I ) 

KAKI I ) = BVI I ) 

KAI I ) = 1 

RETURN 

END 
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SIBFrC BDRY3^ DEBUG 


SUBROUTINE 80RY3^( I # IMf SURF ) 

8DRY34 CORRECTS VALUES CO'IPUTED BY HRB WHEN A HORIZONTAL MESH LINE 
INTERSECTS A BLADE 

COMMON /CALCON/MBIMl ,MBIPl,MBOMl,MBOPl ,MMMI,HM1 ,HT,DTLRtOMLRf 

1 PITCH»CPf EXPON,TWW»CPTIP» TGROG» TB I , TBOr L AMBDA ,TWL 1 1 TMIN, I TMAX , 

2 NIP,[MS(2),BV<2).MVUOO),IVa01),lTVUOO,2),TV(lOO,2)t 

3 DTDMVI 100,2) , BET A V(l 00,2) ,MH( 100,2) ,DTDMH ( 100,2 ), BET AH ( 100,2 ) , 

4 RMH( 100,2 ) ,BEH( 100,2 ) ,RM( 100 ) , BE ( 100 ) , DBOM ( 100) , SAU 100) , 

5 AAA(IOO) 

COMMON /RHOS/RHaHB( 1 00, 2 ) , RHOVB ( 100,2) 

COMMON /HRBAAK/H(4) ,Rt A ) ,B(4),KAK(4) ,KA( A) ,RZ,82,IH(A) 

INTEGER BLOAT,AANOK,ERSOR, STRFN,SLCRO,SURVL, AATEMP,SURF,FIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,MV,MV1 Ml 
IHISURF )*1H( SURF ) 4- 1 
IHS= IH( SURF ) 

HI I)=ABS(MV( IM)-MH( IHS,SURF) ) 

R( I )=RHOHB( IHStSURF) 

B( I )=8EH( IHS,SURF) 

KAK( I ) = 8V( SURF ) 

KA( I ) = 1 

RETURN 

END 


IBFTC SOR DEBUG 

SUBROUTINE SOR 

SUR SOLVES THE SET OF SIMULTANEOUS EQUATIONS FOR THE STREAM FUNCTION 
USING THE METHOD OF SUCCESSIVE OVER-RELAXATION 

COMMON /AUKRHO/ A I 2 SOO , ^ , U ( 2500 ) , K ( 2500 ) , RHO { 2 500 ) 

COMMON /INP/GAM, AR , T I P , RHO I P , WTF L , OMEG A , ORF , BET A I , BET AO , 

1 MBI,M80,MM,NBBI ,NBL, NRSP,MR< 50) ,RMSP( 50 ) , BESP ( 50 ) , 

2 BLOAT, AANDK,ERSOR,STRFN,SLCRD, INTVL tSURVL 

COMMON /CALCON/MBIMl,MBIPl,MBOMl,MBOPl,MMMl,HMl , HT , OTLR , OMLR , 

1 PIT:h,CP,EXPON,TWW,CPTIP,TGROG,TBI,TBO,LAMBOA,TWL, I THIN, UMAX, 

2 NIP, IM$( 2 ) ,BV (2) ,MV( 100) , IVUOl ) , UV( 100,2) |TV( 100, 2) , 

3 DTOMV UOO, 2), BET AVI 100,2) ,MH( 100,2) ,OTOMH( 100,2) ,3ETAH( 100,2) , 

4 RMH( 100,2) ,BEHU00,2 ) ,RMI 100) ,BE ( 100) ,0BDM(100) , SAL I 100) , 

5 AAA(IOO) 

INTEGER BLOAT,AANDK,£RSQR,STRFN,SLCRD, SURVL , AATEMP , SURF , F I RST , 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBOA,LMAX,MH,MLE,MR,MSL ,MSP,MV,MVIM1 
AATEMP = AAIMDK 
IF (0RF,GE.2.) 0RF=0. 

IFIORF.GT.l. ) GO TO 50 
ORF =1. 

ORFUPT= 2* 

40 ORFTEM=ORFOPT 
LMAX = 0. 

50 IF (AATEMP. GT.O) WRUE(6,1010) 

ERROR = 0. 

SOLVE MATRIX EQUATION BY SOR, OR CALCULATE OPTIMUM OVERRE L AXAT I ON 
FACTOR 
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IP - 0 
DO 120 

IPU = I V( IM) 

IPL = 1V( 

IT = irV(IM» 1) 

IF(AArEMP.GT.O) WRITE (6,1020) IM,IT 
DO 120 IP=IPU,IPL 
IPl ^ IP- I 
IP2 = IP+l 

C CORRECT IPl AMO IP2 ALONG PERIODIC BOUNDARIES 
IF( IM.GE.M8I .AND. IM.LE.MBO) GO TO 60 
IF( IT.EQ.ITVl IM, I ) ) IPl = IPI+MBBI 
IF( I T.EQ.ITV( IM, 2) ) IP2 = IP2-MBBI 
60 IT3 = IT 
IT4 = IT 

100 1P3 = IPF( IM-1, I T3) 

IP4 = IPF( IM+1, I TA) 

IF(ORF*GT. 1. ) GO TO 110 
C CALCULATE NEW ESTIMATE FOR LMAX 

UMEW = A (IP, 1 )*U ( IPn+A ( IP ,2 )*imP2) ^^A( IP, 3 J ♦UI IP3 M^A( IP,4)*U( IPA) 
IF (UMEW.LT. l.E-25) U(IP) = 0. 

IF (U( IPI.EQ.O* ) GO TO 115 
RATIU = UMEW/U( I P) 

LMAX= AMAXKRATIO.LMAX) 

U(IP) = UNEW 
GO TO 115 

C CALCULATE NEW ESTIMATE FOR STREAM FUNCTION BY SOR 

110 CHANGE = QRF*(KnP)-U(IP)+A| IP,1)*U(1P1 )+A( IP,2 )+U( IP2 ) 4»A( IP, 3)^ 

1 U( IP3)+A( IP,A )«U( IPA) ) 

ERROR= AMAXl(ERR0ft,A3S(CHANGE)) 

U( IP ) = U( IP )KHAMGE 
115 IF { AATEMP.LE.O) GO TO 120 

WRITE (6,10 30) IT,IP,IPl,IP2,IP3,IPA,(A{IP,n,I = l,A),K(IP) 

120 IT = IT^l 
AATEMP = 0 

IF(UKF.GT.l. ) GO TO 130 

□RFOPT= 2./I l. + SQRM ABS( l.-LMAX) ) ) 

WRirE(6, 1000) URFOPT 

IF(URFIEM-ORFOPT*GT..OOOOl.OR.ORFOPT.GI.1.999) GO TO AO 
WRITE (6,1070) 

ORF = URFOPT 
GO TO 50 

130 IF(ERSOR.GT.O) WR I IE ( 6, 1 OAO ) ERROR 
IF ( ERROR. GT. .000001 ) GO TO 50 
IF(STKFN.LE.O) RETURN 
C 

C PRINT STRrAM FUNCTICN VALUES FOR THIS ITERATION 
C 

WRITE (6,1050) 

DO lAO IM*1,MM 
IPU ^ I V( IMI 
IPL ^ I V( IM+1 )-l 
ITVU = ITV( IM,1) 

WRITE (6,1020) IM,ITVU 
lAO WRITE (6,1060) ( U ( I P ) , I P== I PU , I PL ) 

RETURN 

1000 F0RMAT(2AH ESTIMATED OPTIMUM ORF =,F9.6) 

1010 FORMAT (82HI IT IP IPl IP2 IP3 I PA A(l) A(2) 

1 A(3) A(A) K) 

1020 F0RMAT(5HKIM =, I A , 6X , 64 1 T 1 = ,IA) 

1030 FORMAT! IX, IA,5I6,5F10.5) 
lOAO F0RMAT(8H ERROR =,F11.8) 

1050 FORMAT! IHI,10X,22HSTREAM FUNCTION VALUES) 

1060 FORMAT (2X,10F13.8) 

1070 FORMAT ( IHI) 

END 
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SIBFTC SLAX 


DEBUG 


SUBROUTINE SLAX 

SLAX CALLS SUBROUTINES TO CALCULATE RHO*W-SUB-M THROUGHOUT THE REGION 
AND ON THE BLADE SURFACES, AND TO CALCULAIE AND PLOT THE 
STREAMLINE LOCATIONS 

COMMON / AUKRHO/ A( ^500, A ) , U( 2500 ) , K ( 2500 ) , RHO I 2 500 ) 

COMMON /INP/GAM, AR, T IP, RHOIP , WTFL ,OMEGA,ORF, BETA I ,BETAO, 

1 MBI,MBO,MM,NBBl ,NBL, NRSPfMRI 50) fRMSPI 50) ,BESP(50 J , 

2 BLOAT, AANDK,ERSOR,STRFN,SLCRO, INTVL,SURVL 
COMMON /CALC0N/MBIMI,M3IP1 ,MB0M1 ,MBOPl ,MMM1 ,HMl , HT , OTLK, OMLR , 

1 P{T:h,CP,EXPON,TWW,CPTIP, TGROG, TBI , TBO,LAMBOA,TWL, I TMIN, I TMAX, 

2 NIP, IMS! 2 ) ,BV ( 2) ,MV( 100) , IV( 101 ) , ITV( 100,2) ,TV( 100,2) , 

3 DTOMVI 100,2) , BETA V (100,2 > ,MH( 100,2) ,DTOMH( 100,2) , BET AH ( 100 , 2 ) , 

4 KMH( 100,2) ,BEH( 100,2 ) ,RM( 100 ), BE ( 100) ,OBOM( 100) , SAL 1100), 

5 AAA(IOO) 

COMMON /$LA/TSL(600),UINT(6) 

DIMENSION MSL (100) ,KKK( 14) ,P(4) 

DIMENSION W(2500 ),RWM( 2500), BETA (2500) ,WMB ( I 00 , 2 ) , WTB ( 1 00 , 2 ) , 

1 XDOWN(400) ,YACRDS(400J 

EQUIVALENCE (A(l,l),W(l)),(A(l,2) ,RWM( 1 )),(AI1,3),BETA(1)), 

1 ( A| 1 ,4) ,WMB( l)),{A(201,4),WTB(l)),(A(401,4), XDOWNI I ) ) , 

2 (K( I ) ,YACROS( 1) ) 

INTEGER BLDAT,AANDK,ERSOR, STRFN,SLCRD,SURVL, A ATE MP, SURF, FIRST, 

I UPPER, SI, ST, SRW 

REAL K,KAK,LAMBD A,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
DATA (KKKIJ) ,J=4,14,2)/6*1H+/ 

C 

C CALL SLAVP AND SLAVB8 THROUGHOUT THE REGION 
C 

ITVU= 1TV( 1,1) 

ITVL= ITV( 1,2) 

DO 10 IM=l,MBIMl 
10 CALL SLAVP( IMflTVU, ITVL ) 

00 20 IM=MBI,MBO 
1= 0 

20 CALL SLAVBBI IM) 

90 ITVU = ITV(MBOPl,l» 

ITVL = irV(MB0Pl,2) 

DO 100 IM=MB0P1,MM 
100 CALL SLAVPC IM, ITVU, ITVL) 

C 

C PLOT STREAMLINES 
C 

IF (SLCRD.LE.O) RETURN 
DO 110 IM=1,MM 
ilO MSL( IM) = MV( IM) 

KKK( 1) = 7 
KKK(2) = 6 
KKK( 3) = MM 
P( I ) = 1. 

P( 3) = 0* 

P(4) = 0. 

WRITEIS, 1000) 

CALL PLOTMY(MSL,TSL,KKK,P) 

WRITE(6,1010) 

RETURN 

1000 FORMAT ( 2HPT , 50X , 1 6HSTRE AML 1 NE PLOTS ) 

1010 FORMAT (2HPL,40X,70HSTREAMLINES ARE PLOTTED WITH THETA ACROSS THE 
IPAGE AND M DOWN THE PAGE) 

END 
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$njFTC SLAV 


DEBUG 


SUBROUTINE SLAV 

SLAV CALCULATES RHO*W-SUB-M THROUGHOUT THE REGION AND ON THE BLADE 
SURFACES, AND CALCULATES THE STREAMLINE LOCATIONS 

COMMON SRW, ITER, IEND,LER(2),NER(2) 

COMMON /AUKRHO/ A ( 2SOO, ) , UI 2500 ) , K ( 2500 I , RHOI 2500 J 
COMMON /INP/GAM, AR, T I P, RHOI P , WTFL , OMEGA ,ORF , BET A I , BET AO , 

1 MBI , M80,MM,NBBI , N2L , NRSP , MR I 50) ,RMSP( 50 ) ,BESP(50) , 

2 BLDAT, AANDK,ERSOR,STRFN,SLCRD, INTVL ,SURVL 

COMMON /CALC0N/M6IM1 ,MBIPl,M80Ml,MBOPL,MMMl,HMl ,HT,DTLR,OMLR, 

1 PITCH,CP,EXPON, TWW,CPTIP,TGROG, TB I , T BO, L AMBOA , TWL , 1 TM IN , ITMAX 

2 NIP,IMS(2),BV(2),MV(100),IVIlOl),ITVHOO,2),TVIlOO,2l, 

3 DTDMV(100,2),BETAV(100,2),MH(100,2) ,OTDMH ( 100 ,2 ) ,BETAH( 100,2 ) 

A RMH( 100,2),BEH(100,2>,RM(100),8E(100),OBDM(100), SALUOO), 

5 AAAI 100) 

COMMON /SLA/TSL(600),UINT(6) 

DIMENSION TSPI50 ) ,USPI 50 ) ,DUDT( 50) .TINT (6) 

DIMENSION WI2500) ,RWM12500),8ETA(2500) ,WMB(100,2),WTei 100,2) , 

1 XDUWN(AOO) .YACROSI AOO) 

EQUIVALENCE I AI 1 , 1 ) , W I 1 ) ) , ( A 1 1 , 2 ) , RWM 1 1 ) ) , ( A I 1 , 3 ) ,HET A I 1 ) ) , 

1 I AI 1 ,A) ,WMB( 1 )),IA(201,4) ,WT8I1 ) ) , I A I AO 1 , A ) , XDOWN U ) ) , 

2 IKIl ),YACROS( 1) ) 

INTEGER BLOAT,AANDK,ERSOR, STRFN , SLCRD, SURVL, AATEMP, SURF , F I RST , 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBDA,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 

SLAVP CALCULATES ALONG VERTICAL MESH LINES WHICH DO NOT 
INTERSECT BLADES 

ENTRY SLAVPI IM, I TVU, ITVL ) 

LOC= 0 

NSP= ITVL-ITVU<^2 
IP = IV( IM)-1 
DO iO IT=1,NSP 
IP = IP+1 

TSP(IT) = FLOATI IT^ITVU-1)*HT 
10 USPI IT)= Ul IP) 

USPINSP) = USPID + l. 

IP = I VI IM) 

INTO = INTIUI IP)*5. ) 

IF (U( IP) .GT.O. ) INTU=INTU+1 
DO 20 J=l,5 

UINT (J) = FLOAT! INTU)/5. 

20 INTU = INTU+1 

UINT 16) = UINTI I ) 

GO TO 100 

SLAVBB CALCULATES ALONG VERTICAL MESHLINES WHICH INTERSECT BLADES 
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ENTRY SLAVBB( IM) 

LOC= I 

ITVUPl = ITV< 1) 

ITYLMl = irvUM, 2) 

ITVU = ITVUPl-l 
ITVL = ITVLMl+1 
NSP ^ [TVL-ITVU+1 
TSP( 1) = TV( IM, 1 ) 

TSP(NSP) = TV(IM,2) 

USP( 1) = BVd ) 

USP(NSP) = BV(2) 

IP = IV( IMJ-l 
NSPMl = NSP-1 
IF (2.GT.NSPM1 ) GO TO 70 
DO 60 IT-2, NSPMl 
IP = IP^l 

TSPUT) = FLOAT! !T+ITVU-1)«HT 
60 USP( IT) = U( IP) 

70 DO BO 1=1,6 

8U UINTI I) = FLOAT! I-D/5. 

FOR BOTH SLAVP AND SLAVBB# CALCULATE RHO*W-SUB-M IN THE REGION, AND 
RHO*W AT VERTICAL MESH LINE INTERSECTIONS ON THE BLADE SURFACES 

iOO CALL SPLINE(TSP,USP,NSP,OUOT,AAA) 

IT = LOG 
IPU = IVI IM) 

IPL = IVdM^^D-l 
DO llO IP=IPU,IPL 
IT = ir+1 

ilO RWM(IP) = DUDT( IT)«WTFL/BEdM)/RM( IM) 

120 IF (LOC.EQ.O) GO TO 130 

WMB(IM,l) = DUOT( I )«WTFL/8E(IM) /RM! IM) 

WM8dM,2) = DUOT(NSP)tWTFL/BE(IM)/KM(lM) 

RMDTU2 = (RM( IM)*DTDMV( IM, 1) )«*2 
RMDTL2 = (RM( IM)*DTOMV( IM,2) )^^2 
IF ( HMDTU2-GT. 10000. ) WM6dM,l) = 0. 

IF (KMDTL2.GT.1OC0O. ) WMB(IM,2) = 0. 

WMB(IM,l) = ABS( WMB ( IM, I ) )«SQRT( 1.>RMDTU2) 

WMB(IM,2) = A6S( WMB( IM,2 ) )«SQRT( 1.+RMDTL2) 

130 IF (SLCRD.LE.O) RETURN 
NI = 6 

CALL SPLINKUSP, TSP,NSP,UINT,NI ,TINT,AAA) 

00 140 J=l,6 
L= (J-l)«MM^IM 
140 TSUL)= TINT(J) 

IF (IM.EQd) WRI TE(6,1000) 

WRITE! 6, 1010) MV( IM), (UINTI J),TINT(J),J=l,6) 

RETURN 

1000 FORMAK 1H1/30X,22HSTREAMLINE COORD I NATE S/1 7HL M COORDINATE, 

1 3( 7X, lOHSTREAM FN . , I OX, 5HTHE TA ,4X ) // ) 

1010 FORMAT! IX, 7G18. 7/ (19X, 6G 18.7 ) ) 

END 
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$lliFTC T/.i\G 


DECK 


SUf>><CUTINE TAiNG 

TANG CALCULATES RHO*W-Sun-THETA AND THEN RHO+V. THROUGHCUT THE REGION 
AMO ON THE RLAOE SURFACES, AND CALCULATES THE VELOCITY ANGLE, BETA, 
THROUGHCUT THE REGICN 

COMMON SRW, ITER, I END , LER ( 2 ) , NER ( 2 ) 

CORDON /AUKRHO/ A ( 2 500 , 4 ) , U ( 2500 ) , K ( 2 500 ) , HHO ( 2 500 ) 

CONN'UN / INP/GAH, AK, T IP, RHOIP.WTFL ,OMEGA,ORF,BETA I , RET AO, 

1 Ml! I ,MBO ,MM,NBO I ,NBL ,NRSP,HR ( 50) ,RMSP( 50 ) ,BESP (50 ) , 

2 BLDAT, AANDK, ERSOR ,STRFN,SLCRO, INTVL,SURVL 

CONRCN /CALCON/MRIMI ,MBIP1,MB0M1,MR0PI,,MMMI,HM1 ,HT , OT LR , DMLR , 

1 P ITCH, CP, EXPUN, TWW,CPTI P, TGROG, T6I , TBO , L AMBO A , T WL , I TMIN, I TMAX, 

2 NIP, IMS( 2 ) ,DV ( 2) ,MV( TOO) , IV( 101 ) , ITV( 100,2) ,TV( 100,2) , 

3 LTDMV( 100, 2) , RETAVI 100,2) ,MH( 100,2) ,DTOPH ( 100 ,2 ) , RET AH ( 100 , 2 ) , 

4 RMH( 100, 2 ) ,REH( 100,2 ) ,RM( 100 ) ,BE( 100) ,DCOM( 100) , SAL(IOO), 

5 AAA(IOO) 

DIPENSin-N SPM(10C),USP( 1 00 ) , DUDM ( 1 00 ) 

DIMENSION W( 2500 ) ,RWM( 2 500) ,RETA (2 500) ,WPB( 100,2) ,WTB( 100,2) , 

1 Xn0WN(400) ,YACRUS(400) 

EQUIVALENCE (A(1,1),W(1)),(A( 1 ,2 ) ,RWM( 1)),(A(1,3),PETA(D), 

1 (A(1,4),WMB(1)),(A(201,4),WTB(1)),(A(401,4), xnOV,N( 1 ) ) , 

2 (K( 1) ,YACRQS( 1 ) ) 

INTEGER RLQAT , A A NOK , ERSOR , STRFN , SLCRD , SURVL , A AT EMP, SURF.EIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK , L AMED A , L M A X , MH , ML E , I^R , MSL ,M SP ,M V , MV I P 1 
EXTERNAL HL1,8L2 

PEREOKM CALCUIATIONS ALONG ONE HORIZONTAL LINE AT A TIME 
IT = ITMIN 

10 IF ( IT.GT. ITMAX) RETURN 
SI = 0 

ON THE GIVEN HORIZONTAL MESH LINE, FIND A FIRST POINT IN THE REGION 

IF ( I T .GE .0 . AND. I T.LT.NRB I ) GO TO 60 
IM = MfIMl 
20 IM= IM+1 

IF ( IM.GT.MOQ ) GO TU 200 
SURF = 1 

IE( IT.GE. ITV( IM, D.AND. IT.LT.ITVI IM-1,1) ) GC. TO 70 
IFdM.EO.MBOPl.AND.IT.EO.I TV (MBO, 1 )-l . AND. ITVIMBO, 1 ) - I T V ( MBO , 2 ) 

1 +NBBI.Ea.2) GC TU 70 
SURF = 2 

IF( I T.LE. 1 TV( IM, 2 ) .AND. I T.GT.ITVI IM-1,2) ) GO TO 70 
GO TO 20 

FIRST POINT IS ON BOUNDARY A-H 

60 IMl^ I 
IM = 1 

SPM( 1 ) = MV( 1 ) 

USP (1) = U( IT+1) 
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GU TO ^^0 

FIRST PCI’4T IS ON A BLADE SURFACE 

70 SI = SU’<F 
IMl = 

IM.? = IN 

TH = FLOAT ( 1 T )^HT 
MV IM L = MV( IN 1 ) 

IF ( IM.EQ.MBIPl ) MVIMl = M V I Ml + ( MV (I M2 ) -MV 1 M 1 ) / 1 OOC • 

LFR(2) = 5 

BLCC (VIA ROOT) CALL NO. 5 

IF { S l.FQ. 1 .AND. I r l.NE.MBO)CALL ROOT ( M V I M I , MV ( IMI2 ) , TH , BL 1 , DTLR , 

1 ANS,AAA) 

LER(2) = 6 

BLCC (VIA ROOT) CALL ND . 6 

IF(Sl.E0.2)CALL ROU T ( M V 1 M 1 , K V ( I M2) , TH , BL 2 , D T L R , ANS , A A A ) 

IF(Sl.LO. l.ANO. IMl.EQ.MBO) ANS = MV(MBO) 

SPM( IMi ) = ANS 
USP( IMI )= BV( SI ) 

MOVE ALONG HORIZONTAL MFSH LINE UNTIL MFSH LINE INTERSECTS BOUNDARY 

90 1F( IM.LT.MBI .OP. IM.GT.MBO) GO TO 120 
SURF = 1 

IF ( IT.LT. ITV( !M, SURF ) . AND. 1 T.GF . ITV( IM-1,SURE) ) GO TO 140 
SURF = 2 

IF( IT.GT. ITV( IM, SURD.ANn. IT.LF.ITV( IM-USURF) ) GO TO 140 
120 SPM( IM) = MV( IM) 

IP = IPFI IM, IT) 

USP( IM ) = U( IP ) 

IF ( IM.EC.MM ) GO TO 130 
IM= IM+l 
GO TO 90 

FINAL POINT IS ON BOUNDARY D-F 

130 IMT = MM 
GO rC 130 

FINAL POINT IS ON A BLADE SURFACE 

140 ST = SURF 
IMT= IM 

IMTM1= IMT-1 
TH ^ FLOAT! IT)^hT 
MVIMl = MV( IMTMl ) 

IF ( IMTMl. EQ.MBI ) MVIMl = MV I M 1 + ( M V ( I M2 ) -M V I M 1 ) / 1000 . 

LERI 2) = 7 
C BLCC (VIA ROOT) CALL NO. 7 

IF (ST.EQ. 1 . AND. I MT .NE.MB I )CALL ROC T ( MV I M 1 , M V ( IMT ) , TH,BL1 , 

1 CTLR, ANS, AAA ) 

LERI 2) = 8 

C BLCC (VIA ROOT) CALL NO. 8 

IF ( ST.f 0.2 )CALL ROO T ( M V I M 1 , M V ( I M T ) , TH , B L 2 , 0 TL R , ANS , A AA ) 
IFIST.EO.I.AND. IMT. EQ.MBI) ANS = MV(MBI) 


42 



ooo<^ non nnn 


SPI^( IKT ) = ANS 
USP(ir^T)= BV(ST) 

CALCULATE RHU*W-SUB-THETA AND THfN RHO^W AND BETA IN THE REGION 
150 NSP= IM-IHl+l 

CALL SPL INC ( SPM( OU ) ,USP ( IMl ) ,N5P,0UDM( ll^l),AAA(IM)) 

F1RST=1 

IF ( IKl .NE . I ) F I RST= IM^ 

LAST^ KN 

IF ( l^r .NE.MD) LAST=IMTM1 
IF (FIRST. GT.LAST ) GO TO L ?0 
DU 160 I=FIPSr,LAST 
RWT = -DUOPU I )>!'WTFL/BE( I) 

IP = IPF ( I , IT I 
W(IP) = SORT (RWT^*? + KWM ( IP 
160 BETA(IP) = ATANZ(RWT,RWM( IP) )<=S7. 295779 

CALCULATE RHU=*^W 0^' IHF BLADE SURFACES 

170 IF ( IM I .FO. I ) 90 TO 180 

CALL SEARCH ( SPM ( I M 1 ) » S 1 » I HS ) 

ANS = -OUOM( IM ) <^WTFL/BFH( iHSfSl ) 

WTG(IHS»S1) = ABSIAiMSI^'SORTtl. + U/IRMHllHSfSl )=«=DTORH ( IHS,SU )**2 ) 
180 IF ( IKT.r_C.MM ) GO TO 200 

CALL SI ARCH ( SPM ( I M T ) , S T , I HS ) 

ANS = -DUOM( IMT )<WTFL/BEH( iHSt ST ) 

WTE(IHS,ST) = ABS( ANS SORT! 1. + 1 . / (RMHI I HS , ST ) ♦CTDMH ( IHS,ST ) )«*2 ) 
190 GO IC 20 
200 IT = IT^l 
GO TC: 10 
ENE 


FTC SEARCH DEBUG 

SUBROUTINE SEARCH ( D I ST » SURF 1 1 S > 

SEARCH LOCATES THE POSITION OF A GIVEN VALUE OF M IN THE MH ARRAY 

COMMON /CALCON/MBIMl.MBIPlfMBOMItMBOPLfMMHlvHHl tHT«OTLRt DMLRf 
1 PITCHfCPf EXP0N«TWWtCPTIP«TGROG9TBIfTBQ»LAM80A«TWL»ITMIN,ITMAX9 
Z NIP 9 IMSI2)9BVI2)9MV( 100 ) , IV( 101) • ITV ( L00.2 ) »TV( 100*2 ) 9 

3 0T0HV( 100 92)9 BETAVa0092) tMHC 100t2) vDTOMHC 100*2 I 9BETAH( 100*2 ) . 

4 RMH( 100»2)«BEH{ lOO.Z) *RM( 100 ) *BE ( 100 ) *OBOM ( 1 00 ) 9 SAL(IOO)* 

5 AAA(IOO) 

INTEGER BLOAT 9AANDK«ERS0R9STRFN9SLCRD*SURVL 9 AATEMP* SURF 9 FIRST* 

1 UPPER 9 SI 9 ST 9 SRW 

REAL K,KAK 9 LAMBDAtLMAX,MH,MLE 9 MR 9 MSL 9 MSPtMV 9 MVIMl 
DO 10 I»l*100 

IF (ABSCMHn,SURF)-OIST).GT.DMLR) GO TO 10 
IS = I 
RETURN 
10 CONTINUE 

WRITE (6,1000) OIST 9 SURF 
STOP 

1000 FORMAT C 38HL SEARCH CANNOT FIND M IN THE MH ARRAY/7H DIST :*,G14.6, 
IIOX 96 HSURF >9G14.6) 

END 
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$IBFTC VELOCY DEBUG 


SUBftOUTFf^E VELOCY 
C 

C VELDCY CALLS SUBROUTI'^ES TO CALCULATE DENSITIES AND VELOCITIES 
C THROUGHOUT THE REGION AND ON THE BLADE SURFACESt AND IT PLOTS 
C THE SURFACE VELOCITIES 
C 

COMMON /AUKRHO/ A ( 2 500» 4 > , U I 2500 ) , K { 2500 ) , RHO I 2 500 ) 

COMMON /INP/GAM, AR , T I P , RHO IP » WTFL » OMEGA t ORF , BE T A 1 , 3ETA0 , 

1 M8I,MBO,MM,NB0I ,NBL»NRSPfMR( 50) ,RMSP ( 50 ) , BESP ( 50 > . 

2 BLOAT, AANOKtERSOR* STRFN,SLCRO,INTVL,SURVL 

COMMON /CALC0N/MBIM1,MBIP1,MBGM1,MB0P1 ,MMM1,HM1 ,HT,DTLR,DMLR, 

1 PirCH,CP,EXPON,TWW,CPTIP,TGROG,TBl ,TB0,LAMBDA,TWL,ITM1N, ITMAX, 

2 NIP, IMS(2 ) ,eV(2) ,MV{ 100), IV( ion , ITV( 100,2) ,TV( 100,2) , 

3 DTDMVI 100,2) , BET AVI 100,21 ,MH( 100, 2 > , OTOMH I 100 ,2 ) , BET AH I 100, 2 ) , 

4 RMHC 100,2) ,BEH( 100,2 ) ,RM( 100) ,BEI 100) ,OBDM( 100) , SAL ( 100) , 

5 AAA(IOO) 

DIMENSION KKM18) 

DIMENSION WI2500) ,RWM(2500),BETA(2500) , WM8 ( 100 , 2 ) , WTB I 100, 2 ) , 

I XOOl^NI AOO) ,YACR0SK00) 

EQUIVALENCE iAU,l)fWm), I A ( 1 , 2 ) , RWM I 1) ) , I A ( 1 , 3 ) ,BETA ( 1) ) , 

1 ( AU, A) tWMBI 1)),<A(201,A),WTBU)) , ( A I ^0 1 , A ) , XDOWN ( 1 ) ) , 

2 (K(l ),YACROS( 1) ) 

INTEGER BL OAT, A ANDK,ERSOR,STRFN,SLCRO,SURVL,AATEMP, SURF, FIRST, 

I UPPER, SI, ST, SRW 

REAL K,KAK,LAMB0A,LMAX,MH,MLE,MR,MSL,MSP,MV,MVIM1 
DATA K<K( A)/1H«/ ,KKK I 6) / IHO/ ,KKK I B ) / IH + / ,KKK ( 10 ) / IHX/ 

C 

C CALL VELP, VELBB, AND VELSUR THROUGHOUT THE REGION 
C 

CALL VELPI IfMBlMl) 

CALL VELBB(MBI,M80) 

20 CALL VELPIMBOPl, MM) 

CALL VELSUR 

PREPARE INPUT ARRAYS FOR PLOT OF VELOCITIES 

IFISURVL.LE.O) RETURN 
NP2 = 0 
C TANGENTIAL COMPONENTS 
DO 50 SURF=1,2 
NPl = NP2 
IMSS ^ IMS(SURF) 

IF( IMSS.LT. I ) GO TO 40 
DO 30 IHS=1,IMSS 

IF (ABS(DTOMH(IHS,SURF)«RMH(IHS, SURF)). LT.. 57735) GO TO 30 
NPl = NPl+1 

YACROSINPl) = WT8( IHS,SURF) 

XOOWN(NPl) MH(IHS,SURF) 

30 CONTINUE' 

40 KKK(2*SURF4-I) = NPI-NP2 
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50 NP2 = HPl 

C MERlOrONAL COMPONENTS 
00 BO SUKF=1,2 
NPl = :^P2 

DO 60 IM=MBIP1,M80MI 

IF (A8S(0TDMV(IMtSUKF)«'RM( IM)).GT. 1.7321 ) GD TO 60 
NPl = NPl^l 

YACROS(NPl) = WMB(IM,SURF) 

XDOWN(NPl) = «V( IM) 

60 CONTINJE 

70 KKK(2^SURF+5) = NPI-NP2 
80 NP2 = NPl 

PLOT VELOCITIES 

KKK(l) = 1 
KKK(2) = 4 
P * 5. 

WRITE(6, 1000) 

CALL PLOTMY(XOOWN,YACROS,KKK,P) 

WRITE(6« 1010) 

RETURN 

1000 F0RMAT{2HPT,50X, 24HBLA0E SURFACE VELOCITIES) 

lOlO FORMAT (2HPL,37Xt63HVEL0CITY(W) VS. MERIDIONAL STREAMLINE DISTANCE 
KM) DOWN THE PAGE /2HPL/ 

2 2HPL,50X,50H+ -- BLADE SURFACE 1, BASED ON MERIDIONAL COMPONENT/ 

3 2HPL,50X,50H* - BLADE SURFACE I, BASED ON TANGENTIAL COMPONENT/ 

4 2HPLf 50X,50HX - BLADE SURFACE 2, BASED ON MERIDIONAL COMPONENT/ 

5 2HPL f50X,50H0 - BLADE SURFACE 2, BASED ON TANGENTIAL COMPONENT) 

END 
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$1BFTC VEL 


DEBUG 


SUBROUTINE VEL 

VEL CALCULATES DENSITIES AND VELOCITIES FROM THE PRODUCT OF 
DENSITY TIMES VELOCITY 

COMMON SRW, ITER, IEND»LER(2),NER(2) 

COMMON /AUKRHO/ A( 2 500, 4 ) , U( 2500 ) , K ( 2500 ) , RHD ( 2 500 ) 

COMMON /INP/GAM, AR, TIP,RHOIP,WTFL,OMEGA,ORF,BETAI ,RETAO, 

1 MbI,MBO,MM,NBBI , NBL , NRS P , MR ( 50 ) ,RMSP ( 50 ) , BESP ( 50 ) , 

2 BLDAT,AANDK,ERSOR,STRFN,SLCRO, INTVL fSURVL 

COMMON /CALC0N/MBIM1,MBIP1,MB0M1,MB0P1 ,MMM1,HM1 ,HT,OTLR,OMLR, 

1 PITCH,CP,EXPON, rwWfCPTI P, TGROG, TB I , TBO, LAMBDA ,TWL , ITMIN, I TMAX , 

2 NIP, IMS{2 ) ,BV(2 ) ,MVI100) , tVI 101) , ITVI 100,2 ) ,TV( lOOf 2) , 

3 DTDMVI 100,2) , BETAVI 100,2) ,MH( 100,2) , DTOMH ( 1 00 , 2 ) ,BETAH (100,2 ) , 
^ RMH( 100,2) ,BE HI 100,2 ) ,RM( 100), BE (100) ,D BOM (100), SAL(IOO) , 

5 AAA(IOO) 

COMMON /RHOS/RHOHB( 100,2 ) ,RHOVB ( 100,2 ) 

DIMENSION WWCRMI 100,2) , WWCRT ( 100,2 ) , SURFL( 100,2 ) 

DIMENSION W(2500 ) ,RWM(2500 ), BETA (2500) ,WMB ( 100 , 2 ) , WTB ( 1 00 , 2 ) , 

1 XDOWN(^OO) ,YACR0S(400 ) 

EQUIVALENCE (A(l,l),W(l)),(A(l,2),RWM(l)), ( A ( 1 , 3 ) , BET A { 1)) , 

1 (A(1,^),WMB(D), (A(201,^),WTB(1) ) , ( A ( ^0 1 , ^ , XDOWN ( I ) ) , 

2 (K( 1 ) ,YACROS( 1 ) ) 

VELP CALCULATES ALONG VERTICAL MESH LINES WHICH DO NOT 
INTERSECT BLADES 

ENTRY VELP(FIRST,LAST) 

INTEGER RLOAT,AANOK,ERSOR, STRFN,SLCRD,SURVL, AATEMP,SURF,FIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMBOA,LMAX, MH,MLE,MR,MSL ,MSP,MV,MVIM1 
IF(FIRST.GT.LAST) RETURN 

IF (FIRST. EQ. 1-AND. INTVL.GT.O) WRI TE ( 6 , 1000 ) 

IF (FIRST. EQ.l) RELER = *0 
DO 20 IM*FIRST,LAST 
IPU = I V( IM) 

IPL = IPU+NBBI-l 

TWLMR = 2.»0MEGA«LAMBDA-(0MEGA*RM( IM) )**2 

LERI 1)=4 

DO 10 1P=IPU,IPL 

OENSTY CALL NO. 4 

CALL DENSTY( W( IP ) , RHO ( I P ) , ANS, T WLMR , CP T I P, EXPON , RHO I P , GAM, AR , T I P ) 
10 W( IP) = ANS 

IF ( INTVL.LE.O) GO TO 20 

WRITE (6,1010) IM, ( W( IP ) ,BETA( IP) , IP=IPU, I PL) 

20 CONTINUE 
RETURN 

VEL8B CALCULATES ALONG VERTICAL MESH LINES WHICH INTERSECT BLADES 
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ENTKY VELBB(FIRST,LAST) 

IF(FIftST.GT.LAST ) RETURN 
IF ( FIRST. NE.MBI ) GO TO 30 
SURFLIHBI, 1) = 0. 

SURFL{MBI,2) = 0. 

30 DO 70 IM*FIRST,LAST 
ITVU = ITV(IM,1 ) 

ITVL = ITV( IM,2) 

IPUPl = IPF( IM, I TVU) 

IPLMl = IPF( IM, I TVL ) 

TWLMR = 2.*0MEGA*LAMB0A- (OMEGA»RM( IM) 1**2 
WCR = SQRT(TGROG*TIP*( l.-TWLMR/CPTIP) ) 

IF ( ITVL. LT. ITVU ) GO TO 50 
C ALONG THE LINE BET»<EEN BLADES 
LERI 1 ) = 5 

DO AO IP=IPUPl, IPLMl 
C DENSTY CALL NO. 5 

CALL DENSTYI WI IP ) ,RHO( IP) , ANS, TWLMR, CPTI P, EXPON ,RHO I P, GAM, AR , T I P ) 
40 WI IP) = ANS 

IF I INTVL.LE.O) GO TO 50 

WRITE 16,1010) IM, I WI IP ) ,BETA( IP) , IP=IPUP1 ,IPLM1 ) 

C ON THE UPPER SURFACE 
50 RH08 = RHQVBI IM, 1) 

LERI 1)=6 

C DENSTY CALL NO. 6 

CALL DENSTYIWMBI IM,1),RH0VBI IM, 1 ) ,ANS, TWLMR, CPTIP, EXPON, RHOI P, 

1 GAM, AR, TIP) 

WMBtIM,!) = ANS 
WWCRMIIM,!) = WMOI IM,1)/WCR 
IFI IM.EO.MBI ) GO TO 60 
DELTV = TVI IM-1, 1)-TVIIM,1) 

SURFLIIH,!) = SURFL I IM-1,1 ) + SORTHMVI IM)-MVI IM-1 ) )**2+ 

1 IDELTV*IRMIIM ) +RM I IM- 1 ) ) /2. )**2) 

60 RELER = AMAXl IRELER, ABS I IRH0B-RH0V8I IM, 1 ) ) /RHOVB I IM, 1 ) ) ) 

C ON THE LOWER SURFACE 
RHOB = RHOVBI IM, 2 ) 

LERI 1)=7 

C DENSTY CALL NO. 7 

CALL DENSTYIWMBI I M, 2), RHOVB I IM,2) , ANS, TWLMR , CPT I P,EXPON, RHOI P , 

1 GAM, AR, TIP) 

WMBIIM,2) * ANS 
WWCRMIIM, 2) = WM8I IM,2)/WCR 
IFI IM.EO.MBI ) GO TO 70 
DELTV = TVI IM-1,2)-TVIIM,2) 

SURrHIM,2) = SURFLIIM-1,2)*SQRTI(MV(IM)-MVIIM-1))**2* 

1 IDELTV*IRMIIM)*RMI IM-1) )/2.)**2) 

70 RELER = AMAXIIRELER, ABS I IRH0B-RH0VBIIM,2) I/RHOVBI IM,2) ) ) 

RETURN 

VELSUR CALCULATES ALONG A BLADE SURFACE 

ENTRY VELSUR 
DO 90 SURF=1,2 
IMSS = IMSISURF) 

IF I IMSS.EQ.O) GO TO 90 
DO 80 IHS=1, IMSS 
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TWLMR = 2.*OMEGA+L4MBD4-(OMEGA*RMH( IHSfSURF) )**2 
WCR ^ SQRT(TGROG*riP*(l.-TWLMR/CPTIPn 
RHOB = RHOHBdHSfSURF) 

LER( 11=8 

C DENSTY CALL NO. 8 

CALL DENSTY ( WTB( IHS.SURF )» RHOMB ( I HS , SURF ), ANS , TWLMR ,C PTI P , 

I EXPON,RHOIP,GAM,AR,TIP) 

WTB( IHS, SURF ) = ANS 
WWCRTdHSf SURF) = WTB d HS, SURF ) /WCR 
80 RELER = AMAX1(RELER,ABS( (RH0B-RH0H6( IHS, SURF) )/RHOHB(IHS, SURF) ) ) 

90 CONTINUE 

IF ( RELER.lt. . 00 I ) lEND = lENO+l 
WRITE(6,1080) ITER, RELER 

WRITE ALL BLADE SURFACE VELOCITIES 

IF (SURVL.LE.O) RETURN 
WRITEI6, 1020) 

WRITEI6, lOAO) (MV( I M > , WMB d M , I ) , BETAV (I M, 1) , SURFL ( I M, I ) , 

1 WWCRMdM, l),WMB( IM,2),8ETAV( IM,2),SURFL(1M,2) ,WWCRMdM,2) , 

2 IM=MBI,MBO) 

WRITEI6, 1050) 

DO 100 SURF=1,2 
IMSS = IMS(SURF) 

IFdMSS.LT.l) GO TO 100 
WRITEI6, 1060) SURF 

WRITE I 6, 10 70) (MH( IHS,SURF ) ,WTB ( IHS.SURF ) ,BETAH(IHS,SURF ) , 

I WWCRT ( IHS, SURF ) , IHS=1,1MSS) 

100 CONTINUE 

RETURN 

1000 FORMAT! 1H1/40X,3AHVEL0CI TIES AT INTERIOR MESH POINTS ) 

1010 FORMATdHL,3HIM=,U,5(2AH VELOCITY ANGLEIDEG))/ 

II 5X, 5(G15.4,F9.2 ) ) ) 

1020 FORMAT! 1H1/16X, IH*, 18X,A9HSURFACE VELOCITIES BASED ON MERIDIONAL C 
lOMPONENTS, AOX, IH*/16X, 1H*,53X, 1H*,53X, IH*/ 1 6X , IH* , 19X , 1 5 HBL ADE SUR 
2FACE 1, 19X, 1H*,20X, 15H3LAOE SURFACE 2, 18X, IH*/ 7X , IHM, 8X , 1 H* , 2 ! 3X , 
38HVEL0CITY,3X,23HANGLE!DEG) SURF. LENGTH, 5X, 5HW/ WCR ,6X, 1 H* ) ) 

1090 F0RHAT(!1H ,G13.9,3H * , 2 ( G1 2. 9 , F9.2,2G1 5. 4 , 3H *) ) ) 

1050 FORMAT! IH1/3X,99HSURFACE VELOCITIES BASED ON TANGENTIAL COMPONENTS 
1 ) 

1060 FORMAT! //22X, 15HBLADE SURFACE , 1 1 /7X , IHM, lOX , 8HVEL0C I T Y , 3X , lOHANG 
1LE!DEG),3X,5HW/WCR) 

1070 F0KMAT!1H , 2G 13 . 9, F 9 . 2 , G 15 .9 ) 

1080 FORMAT! 19HLITERATI0N ND . , I 3, 3X , 36HMAX I MUM RELATIVE CHANGE IN DENSI 
ITY =,Gl 1.9) 

END 
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$IBFTC SPLINE DEBUG 


SUBROUTINE SPLINE ( X , Y, N , SLOPE, FM ) 

SPLINE CALCULATES FIRST AND SECOND DERIVATIVES AT SPLINE POINTS 
END CONDITION - SECOND DERIVATIVES ARE THE SAME AT END POINT AND 
ADJACENT POINT 

COMMON 0/80X/S( 1 00 ) , A U 00 ) , B UOO ) , C ( 1 00 ) , F 1 1 00 ) , W( 100) ,S8( 100) , 

1 G(2D0) 

DIMENSION X(N),Y(N) , EM ( N ), SLOPE ( N ) 

INTEGER Q 
DO 10 1=^2, N 
10 SC I )^X( I )-X( I-l ) 

NO=N-i 

IFCNO.lt. 2) GO TO TO 
DO 20 1=2, NO 
AC I )=S( I )/6. 

B( n = CSC I H^SC I^l ) )/3. 

C( I )=S< I + D/6. 

20 F(n = CYCI + l)-Ym)/SIUl)-CYm-YCI-l))/S(I) 

3U ACN) = -.5 
0( 1)=1. 

BCNI^l. 
cm = -.5 
fc L)=0. 

FCN)=0. 

W( l)=B< I ) 

SBC n=cc i)/w( 1 ) 

GC l)=0. 

DO 40 1=2,N 

WC n=8C I )--AC 1 ) + S8C I-l) 

SBC I )=CC I )/W( I) 

40 GCI) = CF( I )-Am*GC I-l > ) /WC I ) 

EM(N)=G(N) 

00 50 1=2, N 
K = N<'I-I 

50 £M(K)=GCK)-SB(K)^EMCK+1 ) 

SLOPEC I )=-S C 2)/6.^C2.»EMC I )+EM( 2) )<^(Y(2)-Ym)/SC2) 

00 60 I=2,N 

60 SLOPEC I )=SC I) /6.«C2.»EMC I ) +EM ( I - U ) + I Y ( I )-Y(f-l))/S(I) 

IF CQ.EQ.13) WRI TE(6,1000) N , C X (I ) , Y ( I) , SLOPE Cl) ,EM C ! ) , I = I , N ) 
RETURN 

1000 FORMAT (2X,15HND. OF POINTS =, I 3/1 OX , IHX, 1 9X , IH Y , 19X , 5HSL0PE , 1 5X , 
12HEM/CAF20.8) ) 

END 
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SUBROUTINE SPLINT ( X , Y , N i Z » M AX , Y I NT , DYOX ) 

C 

C SPLINT CALCULATES INTERPOLATED POINTS AND DERIVATIVES 
C FOR A SPLINE CURVE 

C END CONUITIQN SECOND DERIVATIVES ARE THE SAME AT END POINT AND 
C ADJACENT POINT 
C 

COMMON Q/80X/S( TOO) , A ( 1 00 ) , B ( 100 ) , C (100 ) , F ( 100 ) , W( 100 ) , SB ( 100 ) , 
I G( 100), EMC 100) 

DIMENSION X(N) ,Y(N),Z(MAX) ,YINT(MAX) ,OYQX(MAX) 

INTEGER Q 

IF (MAX.LE.O) RETURN 
I I I = Q 
DO 10 1*2, N 
10 S( I )=X( I )-X( I-l ) 

N0*N~1 

IF(N0.LT.2) GO TO 30 
DO 20 1=2, NO 
A( I )=${ I )/6.0 
8( I ) = ($( I )+S( U1 ) )/3.0 
C( I )=S( 1^1 )/6*0 

20 F(n = (Y(i+i)-Ym)/s(n-i)-(Y(n-Yn-i))/sm 

30 A(N) = -,5 
B( 1)=1.0 
B( N)=1*0 
C( 1) = -.5 
F( 1)=0*0 
F(N)=0.0 
W( 1)=B( 1) 

SBU )=C( 1 ) /W( 1) 

G( i)=0.0 

DO 40 1=2, N 

W( n=8U )~A( I J^SBC I-l) 

SB( I )=cc I) /wn ) 

40 G( n = (Fc I )-A( i)*G( i-in /wn ) 

EM(N)=G(N) 

DO so I=2fN 
K=N+1-I 

50 EM(K)=G{K)-SB(K)*EM(K+1) 

DO 140 1=1, MAX 
K = 2 

IF(Z( I)-X( 1) ) 70,60,90 
60 YINT ( I >=Y( 1 ) 

GO TO 130 

70 IF(Z( I).GE.(m*X( n-.l*X(2) )) GO TO 120 
WRITE (6,1000) Z(I) 

Q = 16 
GO TO 120 
80 K = N 

IF(Z( I) .LE.( 1.1^X(N)-.L*X(N-1) ) ) GO TO 120 
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WRITE (6,1000) Z(l) 

0 = 16 
GO TO 120 

90 IF(Z(I)-X(K)) 120,100,110 
100 YINT( I )=Y(K) 

GO rO 130 
no K = K+1 

lE(K-N) 90,90,80 

120 YINT(I) = EMI K- 1 )*( X ( K) -Z ( 1 ) )**3/6./S(K) +EM(K)* (Zn )-X( K-l )) **3/6. 
l/S(K) + (Y(K)/S(K)-EM(K)’CS(K)/6. )*(Z(I)-X(K-l)) + (Y(K-l)/S(K)-EM(K-l) 
2«S(K)/6. )*(X(K)-Z( I ) ) 

130 DYUX( I )=-EM(K-l ) *(X(K)-Z ( I )) ♦♦2/2.0/S(K)+EM( K)*( X(K-1 )-Z ( I ) ) **2/2- 
l0/S(K)»(Y(K)-Y(K-l))/S(K)-(EM(K)-EM(K-l))*S(K)/6.0 
140 CONTINUE 

MXA = MAX0(N,MAX) 

IF(G.E0.16) WRITE(6,1010) N,MAX,(X(n,ym,Z(I),YINTm, DYDX ( I ) , 

1 1= I, MXA ) 

Q = I 1 1 
RETURN 

1000 FORMAT I S4H SPLINT USEO FOR EXTRAPOLATION- EXTRAPOLATED VALUE = , 
1G14.6) 

1010 FORMAT (2X,21HNQ- OF POINTS GIVEN =,I3,30H, NO. OF INTERPOLATED PO 
TINTS =, 13/ lOX, IHX, 19X, IHY, 16X, I IHX- 1 NTcRPOL. , 9X , I LH Y- 1 NTERPOL . , 
28X, 14Hl)YDX-INTERP0L./(5E20-8n 
END 


51 



o o n o 


SIBFIC SPLN22 DEBUG 


SUBKOUriNE SPLN22 ( X , ¥♦ Y IP » YNP , N* SLOPE , EM) 

SPLN22 CALCULATES FIRST AND SECOND DERIVATIVES AT SPLINE POINTS 
END CONDiriON - DERIVATIVES SPECIFIED AT END POINTS 

COMMON Q/BOX/SaOO) ,M 100) ,B( 100) ,C( 100) »F C 100) , W( 100) , SB ( 100 ) , 

1 G(200) 

DIMENSION X(N),Y(N) , EM ( N )» SLOPE ( N ) 

INTEGER 0 
DO 10 1=2, N 
10 Si I)=X( I )-X( 1-1) 

N0=N-1 

IF(N0.LT.2) GO TO 30 
DO 20 1=2, NO 
A( I )=SI I )/ 6 . 

B( I )=(S( I )+SI I+L ) )/3, 

C( I)=S( I + D/6. 

20 F( I ) = (Y( I + l )-Y( I) )/S( H-1 )-( YU )-Yn-L) )/S( I ) 

30 AIN) = S(N)/6, 

8( l)=S(2)/3. 

B(N) = S<N)/3. 

C( U=S( 2)/6, 

F( l)=(Y(2)-Y( 1) ) /S(2)-YIP 
F(N) = YNP-(Y(N)-Y(N-1) )/S(N) 

W( 1)=B( I) 

SB(1)=C(1)/W(I) 

G( l)=F( I )/W( 1 ) 

DO 40 1=2, N 

W( 1 ) = B( I )-A( D’essai 1-1 ) 

SB( I )=C( I )/W( I) 

40 G( I) = ( F ( I )-A( I )*G( I-I ) ) /we I ) 

EM(N)=G(N) 

DO 50 1=2, N 
K=N+1-I 

50 EM(K)=G(K)-SB(K)«EM(K4-l ) 

SLOPE I l)=-S(2)/6.*(2.^EMe 1)+EM(2) ) + (Y(2)-YU))/SI2) 

DO 60 1 = 2, N 

60 SLOPE! I )=S( I )/6.«{ 2.^EM( I)+EM(I-1))+(Y(I)-YII-1))/S(I) 

IF IQ.E0.18) WRITE!6,1000) N , ( X ( I) , Y < I ) , SL OPE (I ) , EM (I ) , I = 1 , N ) 
RETURN 

1000 FORMAT (2X,15HN0. OF POINTS =, 1 3/ lOX, IHX ,19X, IHY, 19X,5HSL0PE , 15X, 
12HEM/(4F20.8) ) 

END 
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iCBFTC ROOT 


DEBUG 


SUBROUTINE ROOT( A,8, Y,FUNCT, TOLER Y , X , OF X ) 

ROOT FINDS A ROOT FOR (FUNCT MINUS Y) IN THE INTERVAL (A,BI 

COMMON SRW,ITER, lENO t LER ( 2 ) , NER C 2 ) 

INTEGER SRW 

IF (SRW.EQ.2l) WRI TE(6, 1000) A,B,Y,TOLERY 
TOLEKX^ (B-A)/1000. 

AB2= (A^B)/2. 

1 = 0 
X=* A 

10 CALL FUNCT(X,FXf DFX, INF ) 

IF (SRW.EQ.2l) WRITEI6, 1010) 1 » X , F X , OF X » I NF 
IF (ABS( Y-FX).LT.TOLERY) RETURN 
IF ( I.GE.IOOO) GO TO 30 
1 = !♦! 

IF (INF.NE.O .OR. OFX.EQ.O.) GO TO 20 
X= (Y-FX)/OFX^X 

IF (X.GE.A .AND. X.LE.B) GO TO 10 
X A^TOLERX^FLOATI I ) 

IF( I.EO.l) X * B 
GO TO 10 

20 IF (X.LT.AB2) X=X+TOLERX 
IF (X.GE.AB2) X=X-TOLERX 
GO TO 10 

30 WRITE(6,1020) LER(2),A,B,Y 
STOP 

1000 FORMAT I32HIINPUT ARGUMENTS FOR ROOT -- A =G 1 3 . 5 , 3X , 3H3 =»G13.5, 

1 3X,3HY =»G13.5,3X,8HTOLERY =,G13.5/17H ITER. NO. X,17X, 

2 2HFX, 15X*3H0FX, 10X,3HINF) 

1010 FORMAT ( 5X, I3,G16.5,2G18.5, 161 

1020 FORMAT ( lAHLROOT CALL N0..I3/A7H ROOT HAS FAILED TO CONVERGE IN 10 
100 ITERATIONS/AH A = , G1 A . 6 1 1 OX, 3HB =,G1A.6, lOX, 3HY =,GIA.6) 

END 
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SIBFTC BLCD 


DEBUG 


SUBROUTINE BLCD 

BLCD CALCULATES BLADE THETA COORDINATE AS A FUNCTION OF M 
COMMON SRWtITER, I END » LER ( 2 ) , NER ( 2 ) 

COMMON /INP/GAM, AR,TIP,RHOIP,WTFL,OMEGA,ORFtBETAI,BETAO* 

1 MBI»MBO,MM,NBOI ,NBL,NRSPf MR( 50) ♦RMSP( 50) »BESPI50) , 

2 BLOAT, AANDK,ERSOR,STRFN»SLCRDt INTVLfSURVL 

COMMON /CALCON/MBIM1,MB1PI»MBOMI,MBOPI,MMM1,HMI ,HT,DTLR,DMLR, 

1 PITCH,CP,EXPONf TWW,CPTIP, TGROG, TB I , TBO , L AMBD A , TWL , I TMI N , I TMAX , 

2 NIP, IMS(2 } ,BV( 2) ,MV< 100), IV( 101 ) , ITVI 100,2 ) , T Vi 100,2) , 

3 DTDHVi 100,2) , BETA V ( 1 00 , 2 ) ,MH ( 100 , 2 ) , DTDMH ( 100 ,2 ) , BET AH( 100,2) , 

4 RMH( 100,2) ,BEH (100,2 ) , RM I 100 ) , B E ( 100) ,DBOM( 100} , SAL ( 100) , 

5 AAA(IOO) 

COMMON /GEOMIN/ CH0RD(2 ) , STGR( 2 ) , MLE ( 2 ) , THLE i 2 ) , RMl ( 2 ) , RMO( 2 ) , 

1 RI<2) ,R0(2 ),0ETI (2),BET0<2) ,NSPI (2) ,MSP(50,2 ) ,THSP(50,2) 

COMMON /BLCOCM/ EM ( 50, 2 ) , 1 NI T ( 2 } 

ENTRY Btl(M,THETA,OTOM, INF) 

INTEGER BLDAT,AANDR,ERSOR, STRFN,SLCRD,SURVL,AATEMP, SURFtFIRST, 

1 UPPER, SI, ST, SRW 

REAL K,KAK,LAMRDA,LMAX,MH,MLF,MR,MSL,MSP,MV,MVIM1 
REAL M,MMLE,MSPMH,MHM$P 
SURF= 1 
SIGN= 1. 

GO TO 10 

ENTRY BL2(M,THETA,DT0M, INF) 

SURF= 2 
SIGN=-1. 

10 INF= 0 

NSP=* NSPI(SURF) 

IF { INI nSURF ).EQ.13) GO TO 30 
INIT(SURF)= 13 

INITIAL CALCULATION OF FIRST AND LAST SPLINE POINTS ON BLADE 

AA * 0ETI(SURF)/57. 295779 
AA = SIN(AA) 

MSPU,SURF) * RI(SURF)^( l.-SIGN*AA) 

BB = SQRTi l.-AA**2) 

THSPdfSURF) X SIGN*BB*RI(SURF)/RMI (SURF) 

BETKSURF) = AA/ BB/RM I ( SURF ) 

AA X BETO(SURF)/57. 295779 
AA X SIN(AA) 

MSPCNSP,SURF) X CHORD(SURF)-RO(SURF)*( l.+SIGN^AA) 

BB X SQRT( l.-AA»^2) 

THSP(NSP,SURF) X STGR(SURF)^SICN*BB*RO(SURF)/RMO(SURF) 

BETO(SURF) X AA/BB/RMO(SURF) 

DO 20 IAxl,NSP 

MSP( IA,SURF)x MSP{ lA, SURF )+MLE( SURF) 

20 THSP( IA,SURF)x T HSP ( I A, SURF ) ♦THLE ( SURF ) 

CALL SPLN22(MSP( 1 , SURF ), THSP ( 1 , SURF ), BET I ( SURF ) , BE TO ( SURF ) , NSP , 

1 AAA, EMI 1, SURF) ) 

IF(BL0AT«LE*0) GO TO 30 
IF (SURF.EQ.l) WRITE(6, 1000) 

WRITE(6,1010) SURF 

WRITE (6,1020) (MSP( lA, SURF) ,THSP ( I A ,SURF ) , AAA( I A) ,EM( IA,SURF), 

1 IAxI.NSP) 

BLADE COORDINATE CALCULATION 
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30 KK = 2 

IF (M.GT.MSPI If SURF ) ) GO TO 50 

AT LEADING EDGE RADIUS 

MMLE= M-HLE(SURF ) 

IF (MMLE.LT.-DMLR) GO TO 90 
MMLE= AMAXKO.fMMLE J 

THETA= SORT! MMLE*( 2.*RI ( SURF )-MMLE ) )*SIGN 
IF (THETA. EQ.O. I GO TO 40 
RMM= RI (SURF l-MMLE 
DTDM=» RMM/THETA/RHI ( SURF ) 

THETA= THETA/RMI ( SURF )♦ THLE ( SURF ) 

RETURN 
40 INF= I 

OTDM = I.EIO^SIGN 
THETA= THLE(SURF) 

RETURN 

ALONG SPLINE CURVE 

50 IF (M.LE.MSPIKK, SURF n GO TO 60 
IF (KK.GE.NSP) GO TO 70 
KK = KK^l 
GO TO 50 

60 S= MSP(KK,SURF)-MSP(KK-1,SURF) 

EHKM1= EM(KK-1,SURF ) 

EMK= EM(KKfSURF) 

MSPMM= MSP(KK,SURF»-M 
MMMSP* M-MSP(KK-lfSURF) 

THK= THSP(KK,SURF)/S 
THKMl= THSP(KK-l,SURF)/S 

THETA* EMKMl*MSPMM**3/6. /S ♦ EMK*MMMSP**3/6./S 
I MMMSP + (THKHl-EMKMl*S/6. )*MSPMM 
DTDM= -EMKMl*MSPMM**2/2./S ♦ EMK*MMMSP**2/2. /S 
I EMKMl)*S/6. 

RETURN 

AT TRAILING EDGE RADIUS 

70 CMM= CHORD! SURF )+MLE( SURF )-M 
IF (CMM.LT.-OMLR ) GO TO 90 
CMM= AMAXKO.fCMM) 

THETA* SQRT(CMM* (2.*R0( SURF) -CMMI ) *SIGN 
IF (THETA. EQ. 0. ) GO TO 80 
RMM* RD(SURF)-CMM 
OTDM * -RMM/THETA/RHO(SURF ) 

THETA = STGR(SURF)+THETA/RMO(SURF)+THLE(SURFI 
RETURN 
80 INF= I 

DTDM = -1.E104SI6N 
THETA* THLE(SURF)+STGR(SURF) 

RETURN 
C 

C ERROR RETURN 
C 

90 HRITE(6, 1030) LE R ( 2 ) t M , SURF 
STOP 

1000 FORMAT ( IHl , 13X, 33HBLADE DATA AT INPUT SPLINE 
1010 FORMAT! IHLf 17X, 16HBLA0E SURFACE# 14) 

1020 FORMAT ( 7X # IHM # lOX f 5HTHE T A t lOX t lOHOER I VAT I VE # 
1 (4G15.5) ) 

1030 FORMAT ( 14HLBLCD CALL N0.tI3/33H H COORDINATE 
14H M =,G14.6tlOX,6HSURF =,G14.6I 
END 


+ (THK-EMK*S/6. )* 
+ THK-THKM1-(EMK- 


POINTS) 

5XflOH2ND DERIV. / 

IS NOT WITHIN BLADE/ 


55 



n n o r> 


$IBFTC DENSTY DEBUG 


SUBROUTINE OENST YJRHOWtRHO,VEL,TWLMR»CPTIP,EXPON,RHOIP,GAM, AR,T IP) 

DENSTY CALCULATES DENSITY AND VELOCITY FROM THE WEIGHT FLOW PARAMETER 
DENSITY TIMES VELOCITY 

COMMON SRW.ITER, I END, LER ( 2 I , NER ( 2 I 
VEL = RHOW/RHO 
IF I VEL.NE.O. ) GO TO 10 
RHO = RHOIP 
RETURN 

10 TTIP = l.-( VEL**2 + TWLMR) /CPTIP 
IF (TTIP.LT.O. I GO TO 30 
TEMP = TTIP**(EXPON-l. ) 

RHOT = RHOIP*TEMP*TTIP 

RHOWP = -VEL*^2/GAM*RHOIP/AR*TEMP/TIP*RHOT 
IFIRHOWP.LE.O.) GO TO 30 
VELNEW » VEL+IRHOW-RHOT*VEL) /RHOWP 
IFIABSI VELNEW-VEL)/VELNEW.LT..OOOl) GO TO 20 
VEL = VELNEW 
GO TO 10 
20 VEL = VELNEW 
RHO = RHOW/VEL 
RETURN 

30 TGROG = 2.*GAM*AR/(GAM+1. I 

VEL = SORT! TGR0G*TIP*(1.-TWLMR/CPTIP)) 

RHO = RHOIP*( l.-l VEL**2+TWLMR)/CPTIP)**EXP0N 
RWMORW = RHOW/RHO/VEL 
NER( 1) = NER( n + l 

WRITE (6. 1000) LERI 1 ) ,NER I 1 ) .RWMORW 
LFINERI 1 ) .EQ.50) STOP 
RETURN 

1000 FORMAT! 16HLDENSTY CALL N0..I3/9H NERCl) =,I3/10H RHO*W IS ,F7.4, 
134H TIMES THE MAXIMUM VALUE FOR RHO*W) 

END 


ilBFTC IPF DEBUG 

FUNCTION IPFIIM.IT) 

COMMON /CALCON/MBIMl.MB IPl ,MBOM1,MBOP1,MMH1,HM1 ,HT,DTLR,DMLR, 

1 PITCH,CP,EXPON,TWW,CPTIP, TGROG, TB I , TBQ, LAMBDA ,TWL , I TMI N, 1 TMAX , 

2 NIP, IMSI 2),BV(2) ,MV( 100), IVI 101),ITV(100,2) ,TV(100,2), 

3 DTOMVI 100,2) , BET AVI 100, 2) ,MH ( 100, 2 I , DTDMH 1 100,2 ) , BE T AHI 100,2 ) , 

4 RMH(100,2),BEH(100,2),RM(100),BE(100),DBDMI100), SALIlOO), 

5 AAA(IOO) 

IPF = IVI IM)+IT-ITVI IM, I ) 

RETURN 

END 


Lewis Research Center, 

National Aeronautics and Space Administration, 
Cleveland, Ohio, December 9, 1968, 
491-05-00-02-22. 
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